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Abstract 


The  Full  Spectrum  Correlated  ^-distribution  (FSCK)  method,  originally  developed 
for  applications  in  combustion  systems,  is  adapted  for  use  in  shortwave  atmospheric 
radiative  transfer.  By  weighting  /^-distributions  by  the  solar  source  function,  the 
FSCK  method  eliminates  the  requirement  that  the  Planck  function  be  constant 
over  a  spectral  interval.  Consequently,  integration  may  be  carried  out  across  the 
full  spectrum  as  long  as  the  assumption  of  correlation  from  one  atmospheric  level 
to  the  next  remains  valid.  Errors  resulting  from  the  lack  of  correlation  across 
the  full  spectrum  are  removed  by  partitioning  the  spectrum  at  a  wavelength  of 
0.68  /im  into  two  bands.  The  resulting  two-band  approach  in  the  FSCK  formalism 
requires  only  15  quadrature  points  per  atmospheric  layer.  This  represents  a  40- 
90%  reduction  in  computation  time  relative  to  existing  correlated  /^-distribution 
models. 

The  two-band  FSCK  approach  is  developed  for  general  atmospheric  conditions 
through  the  use  of  tabulated  gas  /^-values,  with  nongray  cloud  absorption  coeffi¬ 
cients  added  on-the-fly.  A  two-part  evaluation  of  the  FSCK  calculations  is  pre- 

iii 


sented.  First,  the  two-band  FSCK  results  are  compared  with  line-by-line  (LBL) 
benchmarks  alongside  results  from  an  earlier  radiative  transfer  model  intercompar¬ 
ison  study.  The  median  of  24  1-D  models  included  in  the  intercomparison  has  a 
clear-sky  mean  bias  error  of  -0.27  K/day  relative  to  LBL  benchmark  heating  rates, 
while  the  operational  FSCK  model  has  a  mean  bias  error  of  -1-0.04  K/day.  In  a  sec¬ 
ond  set  of  calculations,  two-band  FSCK  results  are  compared  with  those  from  six 
popular  state-of-the-art  operational  and  research  radiative  transfer  models.  The 
clear-sky  RMS  heating  rate  errors  for  three  empirically-based  models  range  from 
0.78  to  6.28  K/day,  while  RMS  errors  for  three  correlated  /^-distribution  models 
range  from  0.85  to  2.83  K/day.  For  the  same  clear-sky  case  the  FSCK  RMS  error 
is  0.57  K/day.  Cloudy-sky  cases  show  that  the  correlated  /;-distribution  mod¬ 
els  overestimate  in-cloud  heating,  while  the  FSCK  approach  with  nongray  cloud 
absorption  is  closer  to  the  benchmark. 
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Chapter  1 


Introduction 

1.1  A  General  Description  of  Radiative  Transfer 
Models 

The  sun  is  the  external  source  of  energy  that  drives  atmospheric  motions.  As  solar- 
emitted  energy  encounters  the  Earth  atmosphere,  some  of  the  energy  is  absorbed  by 
atmospheric  gases,  some  is  absorbed  by  cloud-liquid,  cloud-ice,  or  aerosol  particles, 
some  is  absorbed  by  the  Earth  surface,  and  the  rest  is  scattered  back  into  space. 

If  we  want  Numerical  Weather  Prediction  (NWP)  models  or  General  Circu¬ 
lation  Models  (GCMs)  to  predict  accurately  a  thunderstorm  this  afternoon,  a 
drought  next  summer,  or  the  average  surface  temperature  100  years  from  now, 
these  models  must  accurately  account  for  the  heating  that  results  from  absorp¬ 
tion  of  solar  energy  by  gases,  clouds,  aerosols,  and  the  surface.  The  models  must 
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also  accurately  represent  the  spatial  and  temporal  distribution  of  the  heating  in 
the  atmosphere  and  on  the  surface.  To  calculate  heating  due  to  absorption  of 
solar  energy  and  the  spatial  and  temporal  distribution  of  that  heating,  NWP  mod¬ 
els  and  GCMs  use  radiative  transfer  parameterizations  which  have  two  essential 
components:  a  solver  and  a  spectroscopic  model. 

A  solver,  as  its  name  implies,  numerically  solves  the  radiative  transfer  eqnation. 
Within  a  NWP  model  or  GCM  the  continnous  vertical  structnre  of  the  atmosphere 
is  typically  approximated  by  discrete  layers,  and  the  horizontal  domain  is  subdi¬ 
vided  into  discrete  colnmns.  There  are  a  variety  of  approaches  to  solving  the 
radiative  transfer  eqnation,  so,  natnrally,  there  are  a  variety  of  solvers  available 
for  use  in  radiative  transfer  parameterizations  with  varying  degrees  of  complexity 
and  accuracy.  In  general,  a  solver  for  solar  energy  entering  the  Earth  atmosphere 
attennates,  layer-by-layer,  the  incoming  solar  energy  in  a  wavelength-dependent 
manner.  Attennation  occurs  as  the  result  of  either  absorption  of  energy  in  a  layer 
or  scattering  of  energy  to  a  different  vertical  layer.  In  a  three-dimensional  solver 
one  can  also  have  have  scattering  that  produces  horizontal  transport  of  radiation 
into,  or  throngh,  neighboring  model  columns.  To  determine  the  amonnt  of  absorp¬ 
tion  and  scattering  in  a  layer,  as  well  as  the  directionality  of  the  scattering,  the 
solver  reqnires  optical  properties  for  the  constitnents  in  each  atmospheric  layer. 

A  spectroscopic  model  provides  optical  properties  for  the  solver  for  each  layer 
based  on  NWP  model-  or  GGM-calculated  pressnre  and  temperatnre,  as  well  as 
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gas,  aerosol-,  and  cloud-particle  concentrations  for  each  layer.  The  solver  requires 
absorption  and  scattering  properties,  including  the  directional  dependence  of  scat¬ 
tering,  for  each  constituent  in  each  atmospheric  layer. 

Gases  in  the  atmosphere  of  Earth  have  absorption  coefficients  that  vary  rapidly 
as  a  function  of  wavelength.  A,  changing  by  several  orders  of  magnitude  across 
the  electromagnetic  spectrum.  As  a  result,  the  most  accurate  way  to  perform 
broadband  radiative  transfer  calculations  is  to  divide  the  full  spectrum  into  nearly 
monochromatic  wavelength  intervals  ((5A  ~  /rm)  and  to  calculate  atmo¬ 

spheric  fluxes  and  heating  rates  for  each  of  these  intervals.  Integrating  the  results 
over  all  intervals  across  the  spectrum  leads  to  broadband  fluxes  and  heating  rates. 
This  approach  is  called  line-by-line  (LBL)  radiative  transfer.  Because  of  the  strong 
variability  of  the  gas  absorption  coefficients  across  the  spectrum,  the  LBL  approach 
requires  millions  of  calculations  to  obtain  a  solution  for  the  full  shortwave,  or  long¬ 
wave,  spectrum,  which  makes  this  method  impractical  for  use  in  NWP  models  or 
GCMs. 

For  over  twenty  years  the  /^-distribution  technique  has  been  used  successfully 
as  a  close  approximation  to  integrated  broadband  LBL  calculations  (e.g..  Goody 
and  Yung  1989;  Lacis  and  Ginas  1991;  Fu  and  Lion  1992).  In  the  /^-distribution 
technique,  gas  absorption  coefficients,  k,  across  a  portion  of  the  spectrum  are 
sorted  in  magnitude  to  create  a  /:-distribution  function.  The  /:-distribution  is  then 
transformed  into  a  cumulative  /^-distribution,  g,  which  is  a  smooth  monotonic 
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function  representing  an  equivalent  nondimensional  wavelength.  Both  the  k-  and 
cumulative  /^-distributions  contain  essentially  the  same  information  for  integrated 
radiative  transfer  calculations  as  the  rapidly  varying  absorption  coefficients  versus 
wavelength  (e.g.,  Lion  2002)  because  the  integration  is  independent  of  the  ordering 
of  the  absorption  coefficients.  However,  instead  of  requiring  hundreds  of  thousands 
of  calculations  to  cover  the  relevant  portion  of  the  spectrum  in  flnx  and  heating 
rate  calculations,  the  cumnlative  A:-distribntion  can  be  integrated  with  only  a  few 
quadrature  points  because  of  its  smoothness. 

The  A:-distribntion  method  is  exact  for  absorbing  non-  and  gray-scattering  ho¬ 
mogeneous  atmospheric  layers  as  long  as  internal  and  external  emission  is  approxi¬ 
mately  constant  over  the  spectral  region  for  which  the  A:-distribntion  is  constrncted. 
To  account  for  the  vertical  inhomogeneity  of  the  atmosphere,  in  the  correlated  k- 
distribution  method  it  is  assnmed  that  the  absorption  coefficients  at  different  pres- 
snres  and  temperatnres  are  correlated  (e.g..  Lion  2002).  Essentially,  this  means 
that  the  absorption  coefficient  in  a  spectral  interval,  SX,  at  one  atmospheric  level 
has  the  same  relative  position  in  the  cnmulative  /;-distribntion  as  the  absorption 
coefficients  in  that  spectral  interval  at  every  other  level,  ensuring  that  the  relative 
sorting  of  absorption  coefficients  from  A-space  to  ^'-space  is  consistent  from  level 
to  level.  Conditions  that  cause  the  relative  sorting  to  change  between  atmospheric 
levels  resnlt  in  a  breakdown  of  correlation  and  introdnce  error  into  the  integration 
in  5'-space. 
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A  limitation  of  the  correlated  /^-distribution  technique  is  that  the  source  func¬ 
tion,  i.e.,  the  Planck  function  for  internal  emission  or  the  solar  source  function  for 
external  emission,  must  be  nearly  constant  over  the  spectral  interval  of  interest 
(e.g..  Goody  and  Yung  1989).  Although  the  Planck  function  varies  smoothly,  it 
covers  a  range  of  several  orders  of  magnitude  across  the  solar  spectrum.  Thus,  in 
broadband  radiative  transfer  calculations,  cumulative  /;-distributions  are  typically 
constructed  for  a  number  of  spectral  bands  (AA  ~  10“^-10°  /xm)  over  which  the 
Planck  function  or  solar  source  function  may  be  assumed  to  be  approximately  con¬ 
stant  (e.g.,  Fu  and  Liou  1992;  Mlawer  and  Clough  1997;  Mlawer  and  Clough  1998; 
Kato  et  al.  1999).  Since  integration  of  each  cumulative  /^-distribution  requires  a  few 
quadrature  points,  integration  of  the  full  spectrum  requires  about  40-150  quadra¬ 
ture  points  per  layer  depending  on  the  particular  model  (e.g.,  Fu  and  Liou  1992; 
Mlawer  and  Clough  1997;  Mlawer  and  Clough  1998;  Kato  et  al.  1999),  which  is 
still  a  considerable  computational  savings  over  LBL  calculations. 

1.2  Motivation  to  Develop  a  New  Shortwave 
Spectroscopic  Model 

Even  considering  the  computational  efficiency  of  the  correlated  /;-distribution  tech¬ 
nique,  40-150  radiation  calculations  is  too  computationally  expensive  for  NWP 
models  and  GCMs.  That  is,  even  as  few  as  40-150  radiation  calculations  in  NWP 
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models  and  GCMs  require  substantial  computational  time  relative  to  all  other 
physical  and  dynamical  calculations.  As  a  result,  radiation  subroutines  in  oper¬ 
ational  models  are  usually  invoked  less  often  than  the  timescale  on  which  cloud 
dynamics  and  microphysics  evolve  in  the  models.  If  the  number  of  calculations 
could  be  further  reduced  without  sacrihcing  accuracy,  the  savings  in  computa¬ 
tional  time  could  be  used  in  both  GCMs  and  NWP  models  to  call  the  radiation 
subroutines  on  a  timescale  that  would  allow  a  more  realistic  interaction  between 
radiation  and  clouds.  In  current  model  practices,  computationally  fast  scaled  k- 
distributions  or  empirical  formulations  are  used  in  the  radiation  parameterization, 
generally  at  the  expense  of  accuracy. 

In  contrast  to  the  correlated  /^-distribution  approach,  in  a  scaled  /^-distribution 
approach  a  cumulative  /^-distribution  is  constructed  for  a  single  reference  level. 
Absorption  coefficients  for  other  levels  are  derived  simply  by  scaling  the  quadra¬ 
ture  values  of  the  reference  cumulative  /^-distribution  as  a  function  of  pressure, 
temperature,  and  absorber  amount.  In  addition  to  scaled  /^-distribution  schemes 
there  are  also  a  variety  of  empirical  methods,  ranging  from  gray  absorption  pa¬ 
rameters  tuned  to  narrow-band  LBL  calculations,  to  empirical  methods  that  create 
functions  that  behave  like  /^-distributions  through  exponential  sum  htting.  Each 
of  these  methods  are  intended  to  provide  accuracy,  while  reducing  computational 
cost  as  much  as  possible. 

The  problem  of  computational  efficiency  of  radiative  transfer  calculations  is 
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not  restricted  to  atmospheric  applications.  The  Full  Spectrum  Correlated  k- 
Distribution  (FSCK)  method  was  developed  to  address  the  problem  of  radiative 
heat  transfer  for  combustion  systems  (Modest  and  Zhang  2002;  Modest  2003). 
The  FSCK  method  eliminates  the  need  for  a  constant  source  function  across 
each  spectral  region  to  be  considered  in  a  radiative  transfer  calculation.  In  the 
FSCK  method  the  cumulative  ^-distribution,  g,  is  no  longer  equivalent  to  a  non- 
dimensional  wavelength.  Instead,  by  weighting  the  /^-distribution  by  the  source 
function,  g  becomes  equivalent  to  a  fractional  source  function.  This  modihcation 
permits  sorting  of  /;-coefhcients  across  bands  as  large  as  the  full  spectrum  as  long 
as  the  condition  of  correlation  between  vertical  atmospheric  levels  is  still  met.  By 
eliminating  the  necessity  for  multiple  spectral  bands,  the  total  number  of  calcula¬ 
tions  can  be  reduced  substantially  without  losing  signihcant  accuracy  relative  to 
broadband  LBL  calculations. 

In  the  case  of  a  purely  external  source  of  radiation,  such  as  solar  radiation 
interacting  with  the  atmosphere  of  Earth,  the  method  is  simplihed  because  the 
only  emission  is  at  the  top  boundary.  However,  the  sun  does  not  radiate  precisely 
as  a  blackbody;  that  is,  the  radiation  from  the  sun  interacting  with  the  Earth 
atmosphere  is  not  precisely  represented  by  a  Planck  function  evaluated  at  the 
effective  surface  temperature  of  the  sun.  Since  estimates  of  solar  irradiance  are 
available  (Kurucz  1992),  the  solar  source  function  can  be  used  as  a  weighting 
function  for  the  /:-distribution  rather  than  the  Planck  function.  Because  traditional 
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/^-distribution  models  assume  a  constant  Planck  function  over  each  spectral  interval 
to  be  integrated,  by  explicitly  accounting  for  the  spectral  variability  of  the  solar 
source  function  the  FSCK  method  has  the  potential  to  improve  the  accuracy  of 
calculations  as  well. 


1.3  Outline  of  this  Work 

In  Chapter  2  we  develop  the  FSCK  method  for  atmospheric  applications  and  sim¬ 
plify  the  equations  for  the  solar  radiation  case.  We  present  a  series  of  standard 
nonhomogeneous  atmospheres,  which  are  the  basis  for  the  initial  tests  of  the  FSCK 
algorithm,  in  section  2.2.  The  first  cases  presented  include  only  the  effects  of 
gaseous  absorption  with  zero  surface  albedo.  Later,  we  present  cases  with  nonzero 
surface  albedo,  molecular  ( “Rayleigh” )  scattering  and  low-  and  high-altitude  clouds 
in  three  different  atmospheres,  including  two  different  treatments  of  cloud  absorp¬ 
tion,  i.e.,  nongray  and  gray.  In  each  of  the  cases,  we  compare  FSCK  broadband 
shortwave  fluxes  and  heating  rates  to  LBL  calculations.  Chapter  2  also  includes 
a  discussion  of  our  preliminary  conclusions  regarding  the  suitability  of  the  FSCK 
method  for  further  development  as  an  atmospheric  radiative  transfer  parameteri¬ 
zation. 

In  Chapter  3  we  advance  the  FSCK  method  from  a  model  that  can  be  applied 
to  an  atmosphere  with  pre-specihed  pressures,  temperatures,  mixing  ratios,  and 
cloud  properties,  as  in  Chapter  2,  to  an  operational  model  that  calculates  fluxes 


and  heating  rates  given  any  input  atmospheric  conditions.  We  discuss  some  of  the 
considerations  affecting  our  choice  of  pressures,  temperatures,  and  mixing  ratios 
used  to  calculate  tabulated  absorption  coefficients,  as  well  as  our  choice  of  major 
and  minor  gas  constituents.  Later,  we  describe  our  method  of  interpolating  from 
tabulated  absorption  coefficients  to  ^-values  appropriate  for  the  atmosphere  of 
interest.  In  sections  3. 2-3. 4  we  describe  and  test  several  different  approaches  for 
combining  constituents,  either  during  the  model  development  process  or  on-the-fly 
during  model  runtime.  We  also  assess  the  relative  accuracy  of  treating  absorption 
by  cloud  particles  as  either  gray  or  nongray  in  the  context  of  the  two-band  FSCK 
approach. 

In  Chapter  4  we  evaluate  the  performance  of  the  two-band  FSCK  model  relative 
both  to  LBL  calculations,  and  to  the  performance  of  other  operational  and  research 
radiative  transfer  parameterizations.  The  evaluation  is  carried  out  in  two  stages. 
First,  in  section  4.1  we  compare  the  FSCK  calculations  to  the  calculations  of  24 
one-dimensional  radiative  transfer  models  that  were  included  in  a  recent  model 
intercomparison  (Barker  et  al.  2003).  Second,  we  compare  FSCK  calculations  to 
those  of  six  specihc  state-of-the-art  operational  and  research  models.  In  section 
4.2  we  describe,  in  detail,  the  methods  used  by  each  of  the  six  models  to  generate 
optical  properties.  In  this  evaluation,  we  embed  the  optical  properties  from  the 
FSCK  model,  the  LBL  model,  and  most  of  the  seven  operational  and  research 
models  into  the  two-stream  solver  of  Zdunkowski  et  al.  (1982),  thus  permitting  us 
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to  make  a  direct  evaluation  of  heating  rate  differences  due  to  treatments  of  optical 
properties.  Chapter  4  concludes  with  a  presentation  of  the  results  of  the  evaluation 
of  the  FSCK  model  relative  to  the  LBL  and  seven  operational  and  research  model 
calculations. 

We  summarize  the  results  in  Chapter  5  and  discuss  a  number  of  suggestions 
to  develop  further  the  FSCK  method  for  applications  in  atmospheric  radiative 
transfer. 
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Chapter  2 


The  FSCK  Algorithm 

2.1  Model  Background 

The  correlated  ^-distribution  and  FSCK  methods  are  similar.  In  both  methods 
the  strongly  variable  absorption  coefficients,  k,  are  sorted  into  smooth,  monoton- 
ically  increasing  cumulative  /^-distributions,  which  can  be  integrated  with  a  few 
calculations.  In  the  correlated  /^-distribution  method,  however,  the  spectral  range 
of  the  absorption  coefficients  to  be  sorted  is  limited  by  the  fact  that  the  emission 
source  function  must  be  constant  over  that  range,  while  in  the  FSCK  method  there 
is  no  restriction  on  the  spectral  range  to  be  sorted.  To  demonstrate  the  difference 
between  the  correlated  /^-distribution  and  FSCK  approaches,  derivations  of  the 
two  methods  by  Modest  (2003),  with  slight  changes  in  the  notation  to  make  it 
more  consistent  with  atmospheric  radiation  nomenclature,  are  reproduced  here  for 
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convenience. 


2.1.1  Correlated  ^-distributions 


The  differential  form  of  the  radiative  transfer  equation  for  an  absorbing,  emitting, 
and  scattering  medium  is,  after  Modest  (2003), 


dl{s,n,\) 

ds 


/cs(^(s).  A) 
4:71 


+kMs),\)B{T{s),\) 
I{s,  kl' ,  X)p{(p{s),  kl' ,  kl ,  X)dkl' , 


'  Atv 


(2,1) 


where  s  is  position,  /  is  monochromatic  radiance,  B  is  the  Planck  function,  T  is 
temperature,  k^  and  kg  are  absorption  and  scattering  coefficients,  p  is  the  scattering 
phase  function,  p  is  an  array  of  state  variables,  including  pressure,  temperature, 
and  mixing  ratio,  and  ^  are  solid  angles,  and  A  is  wavelength.  Assuming  that 
the  absorption  coefficients  are  correlated  across  states  99,  the  set  of  wavelengths  A 
for  which  /ca(v^o,  A)  equals  some  value  k,  where  pq  is  an  arbitrary  reference  state, 
must  also  be  the  set  of  wavelengths  A  for  which  ka{<p,X)  equals  some  value  k*, 
where  ip  is  any  arbitrary  set  of  state  variables.  Since  the  mapping  to  k*  depends 
upon  both  p  and  A,  or  alternatively  k,  we  write  k*  as  a  function  kl  =  kl{p,  k). 

This  mapping  is  represented  mathematically  by  multiplying  equation  (2.1)  by 
the  Dirac  delta  function,  5{k  —  ka{po,  A)),  and  integrating  over  a  spectral  interval 
A  A  sufficiently  narrow  that  the  Planck  function,  scattering  coefficient,  and  scat- 
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tering  phase  function  are  approximately  constant  across  the  interval.  The  result 


is 

dlk(s,  hi) 

- +fc(2(s))]4(s,y) 

+KMs),k)f(^,k)B(T{s))  (2.2) 

ks((p(s))  r  ,  ,  , 

+  7  /  h{sM)p{ip{s)MM)dQ!, 

Jin  ^  ^  ^ 

where 

f{^,k)  =  ^j  5{k  -  k^{(£o,\))d\  (2.3) 

is  the  /^-distribution  function  evaluated  at  the  reference  state  (po, 

Ik{s,Q)  =  ^[  I{sM,X)6{k  -  ks,{(po,X))dX  =  I{sM,k)f{(£o,k)  (2.4) 

Jax 

is  the  probability  of  occurrence-weighted  radiance  I{s,Q^,k)  across  the  interval 
AA,  and  the  wavelength  dependencies  for  the  scattering  coefficient,  phase  function 
and  Planck  function  have  been  dropped  as  they  are  constant  across  AA.  The 
cumulative  /(-distribution  function  evaluated  at  the  reference  state,  i.e., 

pk 

9i^:k)=  /  f{<£o,k')dk',  (2.5) 

Jo 

is  equivalent  to  a  non-dimensional  wavelength.  Dividing  Eq.  (2.2)  by  /(^^o,  k) 
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transforms  the  variable  of  integration  from  k  to  g,  and  we  obtain 


diJs,  n) 

- =  -[K{^{s),g)  +  kMs))]Ig{s,Q) 

+kl{^{s),g)B{T{s))  (2.6) 

47r  ^4^  ^  ^  ^  ^ 

where  Ig  =  Ik/f{(po,k)  .  This  operation  is  possible  since  f{ipo,k)  has  no  depen¬ 
dence  on  location  s.  The  total  radiance  over  AA  is  now  the  integral  of  Eq.  (2.6) 
over  g  from  0  to  1.  Because  the  kl{(p{s),g)  vary  smoothly  with  g,  results  from  Eq. 
(2.6)  can  be  integrated  with  very  few  quadrature  points. 

2.1.2  Full- Spectrum  Correlated  ^-distributions  for 
Longwave  Radiative  Transfer 

The  FSCK  method  extends  the  correlated  ^-distribution  approach  by  relaxing  the 
constraint  that  B{T{s),  A)  be  constant  across  the  interval  AA.  As  a  result,  in  the 
FSCK  approach  AA  can  now  be  arbitrarily  large,  even  the  full  spectrum,  as  long  as 
the  scattering  coefficients  and  phase  functions  are  reasonably  treated  as  constants 
across  AA.  This  is  a  significant  enhancement  to  the  correlated  ^-distribution  ap¬ 
proach,  especially  for  spectral  regions  in  which  the  scattering  coefficients  do  not 
vary  strongly  with  wavelength  (e.g..  Fig.  2.1a),  while  cloud  absorption  coefficients 
do  (e.g..  Fig.  2.1b). 

Although  the  focus  of  this  work  is  an  algorithm  for  shortwave  radiative  transfer. 
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Figure  2.1.  Absorption  and  scattering  coefficients  for  three  different  effective  radii,  rgfr, 
of  a  lognormal  distribution  of  liquid  water  drops  with  a  number  concentration  of  300 
cm“^.  (a)  Scattering  coefficients,  (b)  Absorption  coefficients. 
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we  begin  by  developing  the  case  for  longwave  radiative  transfer,  since  the  shortwave 
case  is  essentially  a  simplihcation  of  the  longwave  case.  The  starting  point  for 
the  FSCK  model  for  longwave  radiative  transfer,  as  in  the  case  of  combustion 
applications  or  atmospheric  emission,  is  again  the  transformation  of  Eq.  (2.1) 
using  the  Dirac  delta  function,  S{k  —  ka{<fo,  A)),  but  the  emission  term  is  treated  in 
a  manner  that  allows  it  to  vary  spectrally.  Multiplying  Eq.  (2.1)  by  S{k  —  ka{(po,  A)) 
and  integrating  over  AA,  we  now  have 

dlk(s,  D) 

- =  -KMs).*:)  +*^!(v(s))]4(s,y) 

+k’Ms),  l-)/(r(s),  k)BUT{s))  (2.7) 

ks((p(s))  f  ,  ,  , 

+  7  /  Ik{s,^)p{<f{s),^,Q)d^, 

where  Bj\x{T{s))  is  the  Planck  function  integrated  over  the  spectral  interval  A  A, 

4(s,  D)=  [  I{sM,X)6{k  -  ks,{(po,X))dX,  (2.8) 

J  AX 

with  no  advantage  of  dehning  4(s,  ^)  in  terms  of  f{<po,  k),  and 

/(T(s),  v>o,  k)  =  )  [  B(T(s),  \)Hk  -  \))d\.  (2.9) 

Note  that  to  eliminate  the  requirement  that  the  Planck  function  must  be  constant 
over  the  spectral  interval  being  sorted,  the  /^-distribution  is  redehned  using  the 
Planck  function  as  a  weighting  function;  that  is,  the  /((-distribution  has  gained 
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temperature  dependence  through  the  Planck  function. 

In  terms  of  the  temperature-dependent  /^-distribution  the  cumulative 
/^-distribution  function  becomes 

g{T{s),^,k)  =  [  f{T{s),(£o,k')dk' 

Jo 

f{T{s),^,k*')dk*'  (2-10) 

=  9{T{s),^,k*). 

The  temperature  dependence  of  the  /;-distribution  in  Eq.  (2.9)  is  problematic, 
as  it  introduces  a  dependence  on  location  s  into  the  /^-distribution.  If  we  were 
to  transform  from  A;  to  g'  in  Eq.  (2.7),  we  would  no  longer  be  able  to  dehne 
Ig  =  h/ f(T{s),i^o,k)  as  the  dependence  of  /  on  s  precludes  us  from  bringing  / 
inside  the  derivative  on  the  left  hand  side  of  the  equation.  To  circumvent  this 
problem  Modest  and  Zhang  (2002)  use  an  arbitrary  reference  temperature  T^f  to 
create  a  reference  cumulative  /^-distribution  function: 

pk 

g{Trei,  >£o,k)=  /  /(Tref,  k')dk' 

Jo 

=  r  /(T,ef,^,r')dr'  (2.11) 

Jo 

=  9iT,eu  V,k*), 

where  T^ei  refers  to  the  temperature  at  which  the  Planck  function  is  evaluated  in 
Eq.  (2.9). 
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It  must  be  noted  that  there  can  be  many  sources  of  emission  (e.g.,  B{T{s)), 
B(T{0)),  /sun)  that  can  be  used  to  dehne  the  emission- weighted  /^-distribution. 
However,  only  one  source  may  be  used  to  dehne  Ig  and  in  a  given  problem. 
The  case  of  a  solar  source  of  emission  will  be  discussed  in  Section  2.1.3,  but  here 
we  continue  the  development  for  the  case  of  longwave  emission. 

Dividing  Eq.  (2.7)  by  Eq.  (2.9)  with  the  temperature  T  evaluated  at  Tref,  we 
obtain 

dig  f(<S,  D ) 

- ^  =  -[^I(^(s),f/ref)  +  kMs))]Ig^^,{s,  D) 

+A::(^(s),(7ref)a(T(s),T,ef,(7ref)i?AA(T(s))  (2.12) 

ks{(p{s))  r  ,  ,  , 

An  ^  ^  ^ 

where 

4ref('5>I^)  =  4(s,  fl)//(Tref,^, /c)  (2.13) 

and 

a(T(s),Tref,fi'ref)  =  f{T{s),(po,  k)  /  f  {T,ei,  (fo,  k) 

(2.14) 

=  dg{T{s),(po,  k)/dg{T,e{,^:  k). 

In  Eq.  (2.14)  a(T(s),  Tref,  (^ref)  is  a  nongray  parameter  (Modest  2003)  that  adjusts 
each  local  cumulative  /^-distribution  in  the  emission  term  to  account  for  differences 
between  the  Planck  function  evaluated  at  the  local  temperature  and  the  Planck 
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function  evaluated  at  the  reference  temperature.  For  notational  convenience  we 
replace  g-^ei  with  g  and  write  the  absorption  coefficient  as  kl(T,-ei,(p,  g)  to  indicate 
explicitly  that  the  absorption  coefficients  are  evaluated  according  to  the  local  state 
variables,  ip,  and  are  subsequently  weighted  in  the  ^-distribution  by  a  Planck 
function  evaluated  at  the  reference  temperatnre  Tref.  With  these  changes  we  have 

di„{s,  n ) 

- =  -[K{Tref,^{s),  g)  +  ks{<£{s))]Ig{s,  Q) 

+k:{T,eu  ip{s),g)a{T{s),T,,,,  g)B^,{T{s))  (2.15) 

ksi<p(s))  /■  ,  ,  , 

+  7  /  Ig{s,Q')p{<p{s),Q\Q)dQk 

Att  ^  ^ 

The  npper  bonndary  condition  for  the  longwave  case  is  /^(stop,  A)  =  0  since 
there  is  no  source  of  longwave  emission  above  Stop-  The  lower  bonndary  condition 
is  given  by 

l\0,Q,X)=e{0,X)B{T{0),X)+^—^^^^[  /^(O,  fi',  A)(cos0)dfi',  (2.16) 
^  J2n  ^  ^ 

where  e(0)  is  assnmed  to  be  a  Lambertian,  i.e.,  isotropic,  emissivity  and  T(0)  is 
the  temperature  of  the  emitting  snrface  (i.e.,  location  0). 

Because  the  upper  boundary  condition  is  trivial  we  do  not  develop  it  fnrther. 
However,  the  lower  bonndary  condition  must  be  transformed  to  use  g  as  the  vari¬ 
able  of  integration.  As  with  the  radiative  transfer  eqnation  we  begin  by  applying 
the  assnmption  of  correlation;  that  is,  we  multiply  Eq.  (2.16)  by  S{k  —  ka{p>o,  A)) 
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and  integrate  over  AA  to  obtain 


4(O,fi)  =  e(O)i?AA(T(O))/(T(O),0o,A;)  +  ^^  /  4(0,  ^)(cos0)d^,  (2.17) 

~T^J2-k 

where  -Baa(7"(0))  is  the  Planck  function  evaluated  at  the  surface  temperature  and 
integrated  over  AA,  and  e(0)  is  the  mean  emissivity  over  the  spectral  interval. 

Since  the  surface  of  the  Earth  emits  radiation  at  temperatures  similar  to  those 
of  the  atmosphere,  it  is  appropriate  to  use  the  same  dehnition  for  the  /^-distribution 
as  was  used  for  internal  atmospheric  emission  in  Eq.  (2.9).  Again  evaluating  the 
Planck  function  in  Eq.  (2.9)  at  Ty-^i  and  using  the  dehnition  of  from 

Eq.  (2.11)  we  divide  Eq.  (2.17)  by  /(Tref,0o;^)  to  arrive  at 

4(0,fi)=e(0)i?AA(T(0))a(T(0),T,ef,A;)  +  ^^^  f  /t(0,  ^)(cos0)d^,  (2.18) 
where  we  have  substituted  g  for  g,-ei  for  notational  convenience. 
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2.1.3  Full- Spectrum  Correlated  ^-distributions 
for  Shortwave  Radiative  Transfer 


In  the  special  case  of  no  internal  source  of  emission  the  second  term  of  Eq.  (2.1) 
is  zero,  so  Eq.  (2.1)  becomes 


dl{s,n,\) 

ds 


ks{^{s),\) 

Att 


(2.19) 


I{s,  X)p{(p{s),  ,  X)dkl' , 


'  Att 


The  boundary  conditions  for  an  external  solar  source  are 


T(stop,  A)  =  /L(stop,  ,  A)  ^  g™^(A)<5(  -  ^)  (2.20) 

for  the  downward  radiance  at  the  top  of  the  atmosphere  (i.e.,  location  Stop)  where 
g™^(A)  is  the  solar  irradiance  at  the  top  of  the  atmosphere,  i.e.,  the  solar  source 
function,  and  flgoi  is  the  solid  angle  of  incident  solar  radiation.  The  lower  boundary 
condition  for  the  upward  radiance  at  the  bottom  of  the  atmosphere  (i.e.,  location 
0)  is  given  by 

Ji(0,  fl,A)  =  -  /  Ji(O,fl',A)(cos0)dfl'.  (2.21) 

We  have  assumed  a  Lambertian,  i.e.,  isotropically  scattering,  surface  with  albedo 
a.  Multiplying  these  boundary  conditions  and  Eq.  (2.19)  by  S{k  —  ka{(po,  A))  and 
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integrating  over  AA,  the  broadband  spectral  interval  of  interest,  we  obtain 


and 


where 


-^sun,AA('®top)  jj.)/(-^sun;  ^0i 


(2.22) 


/^(O,  ff)  =  -  /  4(O,f2')(cos0)dff', 
TT  Jstt 


(2.23) 


dik{s, 


ds 


=  -[Ki^is):  k)  +  ks{(p{s))]Ik{s,  ) 


ks{<£{s)) 

Att 


Ik{s,kl')p{(p{s),kl' ,  kl)dkl', 


'  47r 


(2.24) 


'AX 


=  n,k)f{I^^^,^o,k), 


(2.25) 


/(4un,^,  k)  =  -J - ^ ^  [  4un(stop,  ,  \)5{k  -  k^{ipo,  \))d\,  (2.26) 

4un,AA(Stop,A4)  JaX 

and  Isun  Axi^iop,  is  the  integral  of  the  incident  solar  radiation  over  AA,  the 
broadband  spectral  interval  of  interest. 

Dividing  Eqs.  (2.22)-(2.24)  by  f{l}^^,(po,k),  which  is  independent  of  location 
s  and  can  be  taken  inside  the  derivative  on  the  left  hand  side  of  Eq.  (2.24),  we 


obtain 


and 


where 


il{o,n)  =  - [  i^(o,n'){cos9)dn',  (2.28) 

TT  Jstt  ^ 


diJs,  fi) 

- =  -[Ki.v{.s),  k)  +  A;s(^(s))]4(s,  ) 

k  {  is)) 

+  '7  /  ig{sM)p{^{s)M,n)d^, 

Jin  ^  ^  ^ 

/H(s,  SI)  =  n)//(/i„,^,A^).  (2.30) 


The  cumulative  ^-distribution  function  evaluated  with  the  solar  spectral  distribu¬ 
tion  of  radiance  is 

9{lLn,^,k)=  [  f{I^^^,ipo,k')dk'.  (2.31) 

Jo 

For  an  absorbing  atmosphere  and  a  gray  scattering  atmosphere  with  or  without 
absorption,  only  one  spectral  interval  for  the  entire  spectrum  is  necessary  if  the 
assumption  of  correlation  remains  valid.  One  source  of  loss  of  correlation  in  gaseous 
absorption  between  two  regions  of  space  is  the  presence  of  gradients  of  temperature 
in  the  medium,  as  in  the  case  of  a  combustion  system.  For  example,  at  extremely 
high  temperatures  some  gas  absorption  lines,  so-called  “hot  lines,”  which  originate 
from  excited  transitions,  are  present  at  spectral  locations  where  there  are  no  lines 
at  room  temperature.  Another  reason  for  a  breakdown  in  correlation  is  vertical 
gradients  of  opposite  sign  for  mixing  ratios  of  different  absorbing  species,  leading 
to  a  loss  of  correlation  in  the  absorption  coefficients  with  height. 
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Loss  of  correlation  in  the  atmosphere  results  from  vertical  gradients  of  ozone 
and  water  vapor  absorption  (Figs.  2.2a-b),  as  well  as  with  absorption  in  clouds  and 
their  surrounding  cloud-free  layers  (Figs.  2.2b-c).  High  in  the  atmosphere  (Fig. 
2.2a),  the  strongest  water  vapor  lines  are  similar  in  magnitude  to  ozone  continuum 
absorption  near  0.24  pm,  so  these  lines  will  be  sorted  into  nearly  the  same  part 
of  the  cumulative  /^-distribution.  Lower  in  the  troposphere  (Fig.  2.2b),  the  water 
vapor  lines  are  much  stronger  relative  to  ozone  continuum  absorption  near  0.24  pm, 
so  these  lines  will  be  sorted  into  a  different  part  of  the  cumulative  /^-distribution  at 
that  level.  As  described  in  the  Introduction,  this  change  of  the  relative  sorting  of 
absorption  lines  between  atmospheric  levels  results  in  a  breakdown  of  correlation 
and  introduces  error  into  the  integration  in  g'-space. 

An  example  of  the  manner  in  which  the  presence  of  clouds  can  violate  the 
assumption  of  correlation  is  shown  in  Figs.  2.2b-c.  In  Fig.  2.2b,  which  shows 
absorption  for  a  clear-sky  layer,  absorption  at  point  1  has  a  similar  magnitude  to 
that  at  point  3,  and  is  smaller  than  at  point  2.  However,  when  nongray  cloud 
absorption  coefficients  (Fig.  2.1b)  are  added  to  the  same  atmospheric  layer  (Fig. 
2.2c),  absorption  at  point  1  becomes  greater  than  both  points  2  and  3.  So,  if 
clear-sky  absorption  in  this  layer  had  been  correlated  with  clear-sky  absorption  in 
neighboring  layers,  the  assumption  of  correlation  would  break  down  in  the  presence 
of  cloud.  Examples  of  these  problems  will  be  discussed  in  more  detail  in  Section 
2.2  along  with  an  adaptation  of  the  FSCK  method  to  circumvent  them. 
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Figure  2.2.  Combined  absorption  coefficients  of  water  vapor,  carbon  dioxide,  ozone, 
nitrous  oxide,  carbon  monoxide,  methane,  and  oxygen  at  (a)  0.179  hPa  and  248.2  K,  (b) 
841.783  hPa  and  267.0  K,  and  (c)  the  same  as  (b)  except  with  absorption  coefficients 
of  lognormally  distributed  cloud  liquid  water  drops  with  an  effective  radius  of  5.0  /um 
and  a  number  concentration  of  300  cm“^  added  to  the  gas  absorption  coefficients,  (a) 
and  (b)  The  dashed  line  is  a  point  of  reference  to  compare  the  strengths  of  water  vapor 
lines  relative  to  ozone  continuum  absorption,  (b)  and  (c)  The  dotted  line  is  a  point  of 
reference  to  compare  the  relative  strengths  of  points  1,2,  and  3. 
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2.1.4  Sorting 


In  practice,  the  method  of  sorting  the  absorption  coefficients  to  construct  the  so¬ 
lar  source  function-weighted  /^-distribution  and  then  transforming  it  to  make  the 
equivalent  fractional  solar  source  function-weighted  cumulative  /^-distribution,  g, 
is  straightforward.  The  range  of  absorption  coefficient  values  is  divided  into  bins. 
To  construct  an  ordinary  /;-distribution,  we  compute  the  fraction,  f{k)Sk,  of  the 
wavelength  domain  occupied  by  absorption  coefficients  in  each  “/c-bin,”  i.e.,  be¬ 
tween  k  and  k  +  6k  (e.g..  Goody  and  Yung  1989).  When  an  absorption  coefficient 
falls  within  a  bin  range,  the  bin  count  is  incremented  by  the  monochromatic  wave¬ 
length  interval,  6X,  at  that  spectral  location.  Note  that  if  the  wavelength  interval, 
SX,  were  constant  across  the  spectrum,  it  would  be  sufficient  to  increment  the  bin 
counter  by  one  for  each  occurrence  within  a  /c-bin.  To  construct  a  solar  source 
function- weighted  /^-distribution  (Fig.  2.3c),  rather  than  incrementing  the  bin 
count  by  SX  at  that  spectral  location,  the  bin  count  is  incremented  by  the  value  of 
the  solar  source  function  (Fig.  2.3a)  times  SX  at  that  spectral  location. 

The  range  of  values  of  gas  absorption  coefficients  encountered  in  the  Mc- 
Clatchey  mid-latitude  winter  atmosphere  (McClatchey  et  ah  1972)  falls  be¬ 
tween  10“^^-10^  m“^.  This  range  is  divided  into  log^Q(A;)  bins,  with  100  bins 
per  order  of  magnitude,  for  a  total  of  2401  bins.  Once  the  solar  source  function- 
weighted  /:-distribution  (Fig.  2.3c)  is  constructed,  the  transformation  to  the  cumu¬ 
lative  /;-distribution,  g  (Fig.  2.3d),  is  accomplished  by  adding  up  the  /j-bins,  i.e.. 
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Figure  2.3.  (a)  Solar  Source  Function  (SSF),  (b)  combined  absorption  coefficients 

of  water  vapor,  carbon  dioxide,  ozone,  nitrous  oxide,  carbon  monoxide,  methane,  and 
oxygen  at  148.949  hPa  and  217.9  K,  (c)  SSF-weighted  /^-distribution,  and  (d)  SSF- 
weighted  cumulative  /^-distribution. 
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Qi  =  Qi-i  +  fiSki  where  i  is  the  A:-bin  index,  fi6ki  is  the  value  of  the  ^-distribution  in 
the  bin,  and  Qi-i  is  the  value  of  the  cumulative  distribution  up  to  /c-bin  i  —  1. 
Gaussian  quadrature  is  typically  used  to  integrate  the  cumulative  /^-distribution. 
The  specihc  method  of  determining  quadrature  weights  and  abscissas  used  in  this 
work  will  be  discussed  in  Section  2.2. 


2.2  Fixed  Atmosphere  Case  Studies 

In  the  case  study  that  follows,  the  two-stream  solver  of  Zdunkowski  et  ah  (1982) 
is  used  for  the  LBL  flux  and  heating  rate  calculations.  These  LBL  fluxes  and 
heating  rates  represent  the  “truth,”  i.e.,  the  baseline  for  our  comparison.  For 
consistency,  we  also  used  the  same  solver  in  the  FSCK  calculations.  Therefore,  in 
the  comparisons  between  the  FSCK  calculations  and  the  LBL  standard,  differences 
due  to  model  solvers  are  not  a  consideration. 

2.2.1  Case  Study  Atmospheres,  Gases,  and  Clouds 

The  vertical  prohles  of  pressure,  temperature,  water  vapor  and  ozone  used  for 
the  calculations  in  this  case  study  are  taken  from  the  McClatchey  et  al.  (1972) 
mid-latitude  winter,  subarctic  winter,  and  tropical  atmospheres  and  are  interpo¬ 
lated  to  29  layers  with  a  total  of  30  layer  boundaries,  where  the  lowest  21  lay¬ 
ers  are  1  km  thick,  the  next  two  layers  are  2  km  thick,  the  next  hve  layers  are 
5  km  thick,  and  the  top  layer  is  20  km  thick.  The  prohles  of  oxygen,  carbon 


dioxide,  carbon  monoxide,  methane,  and  nitrous  oxide  are  adapted  from  a  ME- 
TEO  597A  class  project  at  the  Pennsylvania  State  University.  The  atmospheric 
prohles  for  the  McClatchey  et  ah  (1972)  mid-latitude  winter,  subarctic  winter, 
and  tropical  atmospheres  with  gas  constituent  volume  mixing  ratios  used  in  this 
work  are  shown  in  Figs.  2. 4-2. 6.  The  gas  absorption  coefficients  are  extracted 
from  the  HITRAN-96-JPL  database  (Rothman  et  al.  1998)  using  the  Line-by- 
Line  Radiative  Transfer  Model  (LBLRTM)  version  5.21  (Clough  et  al.  1992; 
Clough  and  lacono  1995).  Absorption  coefficients  for  the  individual  gas  con¬ 
stituents  are  combined  by  LBLRTM,  so  line  overlap  and  continuum  absorption 
are  explicitly  calculated  prior  to  sorting  the  combined  absorption  coefficients. 

In  the  cases  where  clouds  are  included  the  cloud  optical  properties  are  cal¬ 
culated  using  Mie  theory  according  to  Toon  and  Ackerman  (1981)  for  lognormal 
distributions  of  liquid-  or  ice-water  spheres  with  a  logarithmic  distribution  width 
of  0.38.  In  our  calculations,  the  liquid-water  cloud  has  a  visible  optical  depth 
of  approximately  30,  while  the  ice  cloud  is  optically  thin  with  a  visible  optical 
depth  of  approximately  0.25.  The  liquid  cloud-drop  distribution  used  here  has 
an  effective  radius  of  5.0  pm.  The  liquid  cloud  in  these  cases,  which  we  also  call 
the  low-altitude  cloud,  is  homogeneously  distributed  between  1000-2000  m  of  al¬ 
titude.  Our  ice-,  or  high-altitude,  cloud-drop  distribution  has  an  effective  radius 
of  75.0  pm  and  is  homogeneously  distributed  in  the  layer  between  7000-8000  m. 
The  cloud  optical  properties  used  for  the  monochromatic  line-by-line  calculations 
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Figure  2.4.  Gas  mixing  ratio  and  temperature  profiles  for  the  mid-latitude  winter 
atmosphere  of  McClatchey  at  ah  (1972).  (a)  Gas  mixing  ratios  (ppmv).  From  right  to 
left  at  30  km:  oxygen  (thin  solid  line),  carbon  dioxide  (dotted  line),  water  vapor  (thick 
solid  line),  ozone  (thick  dashed  line),  methane  (thin  dashed  line),  nitrous  oxide  (dash-dot 
line),  and  carbon  monoxide  (dash-three  dots  line),  (b)  Temperature  (K). 
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Figure  2.5.  Gas  mixing  ratio  and  temperature  profiles  for  the  subarctic  winter  atmo¬ 
sphere  of  McClatchey  at  al.  (1972).  (a)  Gas  mixing  ratios  (ppmv).  From  right  to  left 
at  30  km:  oxygen  (thin  solid  line),  carbon  dioxide  (dotted  line),  water  vapor  (thick  solid 
line),  ozone  (thick  dashed  line),  methane  (thin  dashed  line),  nitrous  oxide  (dash-dot 
line),  and  carbon  monoxide  (dash-three  dots  line),  (b)  Temperature  (K). 
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Figure  2.6.  Gas  mixing  ratio  and  temperature  profiles  for  the  tropical  atmosphere 
of  McClatchey  at  ah  (1972).  (a)  Gas  mixing  ratios  (ppmv).  From  right  to  left  at  30 
km:  oxygen  (thin  solid  line),  carbon  dioxide  (dotted  line),  water  vapor  (thick  solid  line), 
ozone  (thick  dashed  line),  methane  (thin  dashed  line),  nitrous  oxide  (dash-dot  line),  and 
carbon  monoxide  (dash-three  dots  line),  (b)  Temperature  (K). 
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were  produced  at  a  spectral  resolution  of  0.01  /iui  and  then  interpolated  to  the 
wavenumber  resolution  of  the  gas  absorption  coefficients. 

2.2.2  Clear-Sky  Calculations 

The  hrst  case  is  a  basic  test  of  the  FSCK  method  versus  LBL  calculations.  In  this 
case  there  are  no  clouds,  no  molecular  scattering,  the  surface  albedo  is  zero,  and  the 
mid-latitude  winter  atmospheric  prohle  is  used.  The  full  shortwave  spectrum  (0.24- 
4.60  /im)  is  treated  as  a  single  spectral  band  and  the  gas  absorption  coefficients 
are  sorted  into  a  single  cumulative  /^-distribution  per  layer.  For  this  initial  test 
of  the  FSCK  method,  8192  quadrature  points  are  used  for  the  integration  of  the 
cumulative  /^-distributions  in  order  to  ensure  that  differences  between  the  FSCK 
results  and  the  LBL  calculations  are  not  due  to  poor  quadrature.  This  number 
was  selected  by  doubling  the  number  of  quadrature  points,  starting  with  one,  until 
the  calculated  heating  rate  prohle  ceased  to  change  with  a  subsequent  doubling  of 
points. 

Comparing  the  vertical  prohles  of  downward  huxes  and  heating  rates  for  the 
LBL  calculations  and  the  FSCK  results  (Figs.  2.7a-c),  we  hnd  less  molecular 
absorption  in  the  FSCK  results  with  the  most  signihcant  differences  in  the  lowest 
14  km  of  the  atmosphere.  The  heating  rate  RMS  error  for  the  full  atmospheric 
depth  is  approximately  4.2%,  but  is  nearly  20%  for  the  lowest  14  km.  As  Figs. 
2.2a-b  illustrate,  these  differences  result  from  a  loss  of  correlation  between  the 
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Figure  2.7.  Upward  and  downward  fluxes  (a,  d,  g,  and  j),  heating  rates  for  all  29  layers 
(b,  e,  h,  and  k),  and  heating  rates  for  the  lowest  14  layers  (c,  f,  i,  and  1)  for  (a-c)  the 
full  shortwave  spectrum  (0.24-4.6  ^m)  treated  as  a  single  band,  (d-f)  band  1  (0.24- 
0.68  ;um),  (g-i)  band  2  (0.68-4.6  /um),  and  (j-1)  the  full  shortwave  spectrum  treated  by 
adding  bands  1  and  2.  The  LBL  solution  ,  the  8192  quadrature  point  solution,  and  the  15 
quadrature  point  solution  are  given  by  the  solid,  dashed,  and  dotted  lines,  respectively. 
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/^-distributions  for  different  vertical  levels  as  a  result  of  different  vertical  gradients 
in  the  water  vapor  and  ozone  mixing  ratios.  This  condition  is  always  present  in 
the  atmosphere  of  Earth. 

Since  absorption  by  ozone  is  dominant  in  the  upper  atmosphere,  less  absorbing 
in  the  lower  atmosphere  relative  to  the  other  gases,  and  nearly  the  only  mechanism 
for  absorption  in  the  spectral  range  from  0.24-0.68  /rm,  we  separate  the  shortwave 
spectrum  into  two  bands,  with  band  1  from  0.24-0.68  /xm  and  band  2  from  0.68- 
4.60  /im,  and  test  the  assumption  of  vertical  correlation  within  each  of  these  two 
spectral  bands. 

This  approach  is  similar  to  the  part-spectrum  approach  described  by  Modest 
and  Riazzi  (2004),  in  which  spectral  regions  with  similar  scattering  properties  are 
grouped  together.  In  our  case  the  spectral  separation  is  not  intended  to  group 
together  parts  of  the  spectrum  with  common  scattering  properties.  Rather,  the 
spectral  separation  isolates  the  part  of  the  shortwave  spectrum  dominated  by  ozone 
continuum  absorption  from  the  remainder  of  the  shortwave  spectrum. 

The  resulting  fluxes  and  heating  rates  are  much  more  accurate  (Figs.  2.7d-l). 
For  band  1  the  maximum  errors  in  the  flux  in  the  lowest  layers  of  the  atmosphere 
are  less  than  1%  and  heating  rate  errors  are  also  small,  being  less  than  approxi¬ 
mately  0.005  K/day  (Figs.  2.7d-f).  We  were  able  to  hnd  a  5-point  quadrature  of 
the  FSCK  cumulative  distribution  functions  for  band  1  that  led  to  nearly  identical 
results  with  the  8192-point  quadrature.  That  is,  the  FSCK  treatment  of  gaseous 
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absorption  in  band  1  is  accurate  to  approximately  3%  of  LBL  heating  rates  and  a 
5-point  quadrature  for  this  band  is  possible  without  a  signihcant  loss  of  accuracy. 

The  band  2  comparisons  (Figs.  2.7g-i)  demonstrate  that  the  assumption  of 
vertical  correlation  is  also  valid  for  the  spectral  interval  from  0.68-4.60  /rm,  as 
the  FSCK  8192-point  quadrature  flux  and  heating  rate  differences  from  the  LBL 
results  are  less  than  1%  and  6%,  respectively.  However,  for  band  2,  the  quality 
of  the  FSCK  results  are  highly  dependent  on  the  number  of  quadrature  points. 
Using  a  10-point  quadrature,  which  represents  a  compromise  between  accuracy 
and  computational  efficiency,  we  found  the  heating  rate  errors  in  the  troposphere 
to  be  less  than  13%,  while  the  errors  in  the  stratosphere  reached  66%.  But,  the 
consequence  of  the  large  percentage  error  in  the  stratosphere  amounts  to  only 
a  heating  rate  error  of  1.1  K/day,  which  is  small,  i.e.,  less  than  6%,  when  the 
total,  full  spectrum  heating  rate  is  considered  (Figs.  2.7j-l).  Overall,  the  8192 
quadrature-point  FSCK  method  produced  heating  rates  accurate  to  approximately 
3.8%  in  the  lowest  14  km  and  3.0%  for  the  full  depth  of  the  atmosphere.  Limiting 
the  two-band  quadrature  scheme  to  5  points  in  band  1  and  10  points  in  band  2  led 
to  larger  full  spectrum  heating  rate  errors  of  5.1%  and  7.3%  in  the  lowest  14  km 
and  the  full  atmospheric  depth,  respectively. 

Shifted  Gauss- Chebyschev  II  quadrature  (Wang  and  Modest  2004)  is  used 
to  calculate  the  weights  and  abscissas  for  band  2  (Table  2.1).  However,  there  is 
a  sharp  kink  present  in  many  of  the  cumulative  /^-distributions  of  band  1  near 
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Table  2.1.  Trapezoidal  quadrature  abscissas  and  weights  for  band  1  and  shifted  Gauss- 
Chebyschev  II  quadrature  abscissas  and  weights  for  band  2  before  and  after  transforma¬ 
tion. 


Band  1  Band  2  before  Band  2  after 

transformation  transformation 

Abscissa  Weight  Abscissa  Weight  Abscissa  Weight 


0.453830 

0.929175 

0.962060 

0.981040 

0.994325 


0.907660 

0.043030 

0.022740 

0.015220 

0.011350 


0.000000 

0.156434 

0.309017 

0.453990 

0.587785 

0.707107 

0.809017 

0.891007 

0.951057 

0.987688 


0.076046 

0.160133 

0.144407 

0.144940 

0.122105 

0.116036 

0.087385 

0.076215 

0.043742 

0.028991 


0.000000 

0.225220 

0.425618 

0.596540 

0.735342 

0.841487 

0.916537 

0.964017 

0.989172 

0.998634 


0.113893 

0.220272 

0.179779 

0.160402 

0.117413 

0.094052 

0.057194 

0.037684 

0.014493 

0.004818 


5'=0.95  (not  shown)  that  reduces  the  accuracy  of  Gaussian  quadrature  methods. 
As  a  result,  we  use  trapezoidal  quadrature  for  band  1  with  weights  and  abscissas 
optimally  selected  to  fit  the  LBL  solution. 

The  next  test  case  also  consisted  of  the  same  clear-sky  mid-latitude  winter 
prohle,  but  with  molecular  (i.e.,  “Rayleigh”)  scattering  and  a  Lambertian  surface 
with  gray  albedo  of  0.15  included  in  both  the  FSCK  and  LBL  calculations.  For 
the  LBL  calculations  we  used  the  spectrally-dependent  Rayleigh  scattering  opti¬ 
cal  depths  calculated  according  to  Teillet  (1990),  while  for  the  two-band  FSCK 
calculations  the  Rayleigh  scattering  optical  depths  were  calculated  using  the  same 
technique  but  at  fixed  wavelengths  of  0.47  pm  for  band  1  and  0.90  pm  for  band 
2.  In  the  FSCK  calculations  we  also  used  the  same  5-  and  10-point  quadrature 
values  for  bands  1  and  2  that  we  derived  in  the  gas-absorption  only  case.  Not 


37 


surprisingly,  the  flux  and  heating  rate  errors  that  we  obtained  for  this  case  (Fig. 
2.8a-c)  do  not  differ  in  any  significant  way  from  those  that  we  obtained  for  the 
gas-absorption  only  case.  Moreover,  results  for  clear-sky  calculations  in  subarctic 
winter  and  tropical  atmospheres  are  similar  (not  shown). 

2.2.3  Cloudy  Sky  Calculations 

The  next  two  sets  of  test  cases  examine  the  accuracy  of  FSCK  calculations  in  cloudy 
conditions.  In  the  first  case,  an  FSCK  treatment  of  nongray  cloud  absorption  with 
parameterized  gray  cloud  scattering  is  evaluated,  and  in  the  second  case,  both 
gray  cloud  absorption  and  gray  cloud  scattering  are  parameterized.  In  both  sets 
of  calculations  we  consider  the  effects  of  low  cloud,  high  cloud  and  both  low  and 
high  clouds  combined.  The  two  sets  of  calculations  differ  in  the  method  by  which 
we  account  for  cloud  absorption  in  the  FSCK  calculations. 

In  the  first  set  of  calculations,  where  cloud  absorption  is  nongray  and  cloud 
scattering  is  gray,  we  interpolated  the  O.Ol-pm  resolution  cloud  absorption  coef¬ 
ficients  (Fig.  2.1b)  to  the  spectral  resolution  of  the  gas  coefficients  in  the  two 
cloud  layers  and  then  added  the  interpolated  cloud  absorption  coefficients  to  the 
gas  absorption  coefficients.  The  combined  cloud  and  gas  absorption  coefficients 
were  then  sorted  to  form  the  cumulative  ^-distribution  functions  from  which  the 
final  5-  and  10-point  quadrature  ^-values  for  the  radiative  transfer  were  selected. 
The  advantage  of  this  approach  is  that  spectral  variations  in  the  cloud  absorption 
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Figure  2.8.  Upward  and  downward  fluxes  (a,  d,  g,  and  j),  heating  rates  for  all  29 
layers  (b,  e,  h,  and  k),  and  heating  rates  for  the  lowest  14  layers  (c,  f,  1,  and  1)  for  (a-c) 
clear  sky,  (d-f)  a  liquid  cloud  with  nongray  absorption  and  gray  scattering,  (g-i)  an  ice 
cloud  with  nongray  absorption  and  gray  scattering,  and  (j-1)  both  clouds  with  nongray 
absorption  and  gray  scattering.  The  LBL  solution,  the  8192  quadrature  point  solution, 
and  the  15  quadrature  point  solution  are  given  by  the  solid,  dashed,  and  dotted  lines, 
respectively. 
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coefficients  are  explicitly  taken  into  account.  The  disadvantage  to  this  approach 
is  that  adding  the  cloud  absorption  coefficients  to  the  gas  absorption  coefficients 
leads  to  a  partial  loss  of  correlation  in  the  /^-distributions  between  the  cloudy  layers 
and  their  neighboring  clear  sky  layers,  as  we  illustrated  in  Fig.  2.2b-c. 

The  parameterized  gray  liquid-  and  ice-cloud  asymmetry  parameters,  single¬ 
scattering  albedos,  and  scattering  cross  sections  used  for  this  hrst  set  of  nongray 
cloud  absorption  and  gray  cloud  scattering  calculations  are  shown  in  Table  2.2. 
We  obtained  these  values  by  taking  a  band  average  of  the  hne-spectral-resolution 
Mie  calculations  weighted  by  spectrally-varying  direct  downwelling  transmission 
at  500  hPa.  The  resulting  parameters  were  hne-tuned  by  minimizing  the  flux 
errors  between  the  8192-point  solution  and  the  LBL  calculations.  Note  that  in 
this  case  nongray  absorption  by  cloud  drops  is  explicitly  accounted  for  by  the  k- 
distributions,  whereas  absorption  by  cloud  drops  is  typically  accounted  for  by  gray 
cloud  extinction  cross  section  and  single-scattering  albedo  parameters.  In  order  to 
avoid  accounting  for  absorption  by  cloud  drops  twice,  the  gray  cloud  parameters 
are  adjusted.  First,  the  gray  cloud  parameter  representing  the  extinction  cross 
section  is  given  the  value  of  the  scattering  cross  section;  that  is,  the  extinction 
cross  section  times  the  single-scattering  albedo.  Next,  the  single-scattering  albedo 
is  set  to  1.0.  This  procedure  essentially  subtracts  the  gray  cloud  absorption  cross 
section  from  the  gray  cloud  extinction  cross  section. 

As  Figs.  2.8d-f  illustrate  for  a  low  cloud  and  Figs.  2.8g-i  illustrate  for  a  high 
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Table  2.2.  Asymmetry  parameter,  g,  single-scattering  albedo,  loq,  and  extinction  cross 
section,  aext,  for  a  lognormal  distribution  of  liquid  water  spheres  of  radius  5  ;um  and  ice 
spheres  of  radius  75  /um  for  band  1  (0.24-0.68  ^m)  and  band  2  (0.68-4.60  /um),  nongray 
absorption  and  gray  scattering,  and  gray  absorption  and  gray  scattering. 


Band  1 

Band  2 

9 

(Jo 

^ext 

9 

Uq 

^  ext 

Nongray  absorption  and 
gray  scattering 

Liquid 

Ice 

0.849 

0.888 

1.000 

1.000 

1.10e-h02 

2.33e-h04 

0.825 

0.898 

1.000 

1.000 

9.40e-h01 

2.35e-h04 

Gray  absorption  and 
gray  scattering 

Liquid 

Ice 

0.856 

0.867 

0.999 

0.999 

1.10e-h02 

2.32e-h04 

0.800 

0.898 

0.998 

0.995 

9.30e-h01 

2.34e-h04 

cloud,  maximum  upward  and  downward  flux  errors  were  always  less  than  6%  with 
typical  errors  near  1%,  while  the  heating  rate  errors  were  generally  less  than  4% 
in  the  clear  sky  region  above  the  highest  cloud  layer  and  no  greater  than  8-12% 
in  and  below  the  cloud  layer.  In  both  the  low-  and  high-cloud  cases  the  FSCK 
simulations  with  nongray  cloud  absorption  and  gray  cloud  scattering  produced 
slightly  less  heating  than  the  LBL  calculations.  The  cloud  parameters  that  we 
obtained  for  the  single  cloud-layer  cases  were  used  for  the  combined  low-  and  high- 
cloud  case,  and  results  of  that  case  were  consistent  with  the  single  cloud-layer  cases 
(Figs.  2.8j-l). 

The  likely  source  of  error  in  this  case  of  nongray  cloud  absorption  is  a  breakdown 
in  correlation  between  the  cloud  and  surrounding  clear-sky  layers,  as  described  in 
Section  2.1,  and  correlation  could  be  partially  restored  by  subdividing  band  2 
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into  more  bands.  However,  we  are  attempting  to  balance  accuracy  with  efficiency. 
Each  additional  band  would  require  5-10  additional  quadrature  points,  and  since 
the  magnitude  of  the  maximum  in-cloud  heating  rate  errors  in  this  nongray  cloud 
absorption  case  are  of  similar  magnitude  to  heating  rate  errors  of  other  state-of-the- 
art  one-dimensional  radiative  transfer  models  relative  to  line-by-line  calculations 
(Barker  et  ah  2003),  we  do  not  attempt  to  reduce  further  the  two-band  FSCK 
errors. 

In  the  second  set  of  calculations  we  treated  the  cloud  absorption  coefficients  as 
gray  within  each  band,  in  addition  to  gray  scattering  coefficients  and  asymmetry 
parameters.  The  values  of  the  three  parameterized  gray  cloud  properties,  i.e., 
asymmetry  parameter,  single-scattering  albedo  and  extinction  cross-section,  were 
determined  as  described  above  for  the  gray  cloud  scattering  properties  in  the  hrst 
set  of  calculations,  except  that  the  gray  cloud  extinction  cross  sections  and  single¬ 
scattering  albedos  now  account  for  cloud  absorption.  The  tuning  parameters  for 
gray  cloud  absorption  and  gray  cloud  scattering  are  shown  in  Table  2.2. 

As  the  results  in  Figs.  2.9a-c  and  d-f  illustrate,  values  for  the  three  coefficients 
exist  that  provide  accurate  fluxes,  to  better  than  3%,  and  heating  rates,  to  better 
than  2%,  for  each  cloudy  atmosphere.  These  results  demonstrate  that  relatively 
accurate  parameterization  of  single  cloud-layer  optical  properties  with  gray  values 
for  the  two  bands  is  possible.  Note  that  with  three  tuning  parameters,  i.e.,  asym¬ 
metry  parameter,  single-scattering  albedo,  and  extinction  coefficient,  flux  errors 
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Figure  2.9.  Upward  and  downward  fluxes  (a,  d,  and  g),  heating  rates  for  all  29  layers 
(b,  e,  and  h),  and  heating  rates  for  the  lowest  14  layers  (c,  f,  and  i)  for  (a-c)  a  liquid 
cloud  with  gray  absorption  and  gray  scattering,  (d-f)  an  ice  cloud  with  gray  absorption 
and  gray  scattering,  and  (g-i)  both  clouds  with  gray  absorption  and  gray  scattering.  The 
LBL  solution,  the  8192  quadrature  point  solution,  and  the  15  quadrature  point  solution 
are  given  by  the  solid,  dashed,  and  dotted  lines,  respectively. 
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might  be  reduced  more  with  a  different  combination  of  these  parameters.  If  this 
gray  cloud  absorption  and  gray  cloud  scattering  approach  is  further  developed,  sets 
of  parameters  must  be  determined  for  the  range  of  cloud  particle  size  distributions 
encountered  in  the  numerical  models  in  which  the  FSCK  approach  is  implemented. 

Incorporating  the  coefficients  computed  separately  for  the  low  and  high  clouds 
into  a  combined  low  and  high  cloud  FSCK  calculation  produced  flux  and  in-cloud 
heating  rate  errors  of  1%  and  6-14%,  respectively,  relative  to  the  LBL  results 
(Figs.  2.9g-i).  The  magnitudes  of  the  in-cloud  heating  rate  errors  using  this  second 
approach  are  similar  to  those  resulting  from  explicitly  incorporating  the  nongray 
cloud  absorption  coefficients  into  the  solar  source  function-weighted  cumulative 
/^-distributions,  although  the  magnitudes  and  signs  of  the  errors  are  less  consistent 
from  the  single-layer  to  two-layer  cloud  cases  than  for  the  nongray  cloud  absorption 
calculations. 

As  we  just  demonstrated,  the  two  sets  of  FSCK  calculations  for  the  mid-latitude 
winter  atmosphere  produced  comparable  results,  with  errors  in  the  fluxes  and  heat¬ 
ing  rates  of  similar  magnitudes.  Performing  identical  sets  of  calculations  for  both 
the  subarctic  winter  and  tropical  atmospheres,  with  the  same  quadrature  points 
and  cloud  parameters  derived  for  the  mid-latitude  winter  atmosphere,  we  obtained 
results  consistent  with  those  from  the  mid-latitude  atmosphere  calculations. 
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2.2.4  Case  Study  Observations 


The  flux  and  heating  rate  results  obtained  in  this  case  study  demonstrate  that 
the  full  spectrum  correlated  /^-distribution  (FSCK)  method  is  an  accurate  and 
efficient  alternative  to  line-by-line  (LBL)  calculations  for  shortwave  atmospheric 
radiative  transfer.  By  weighting  /^-distributions  by  the  Planck  function,  or  by  the 
solar  source  function  in  the  case  of  solar  radiation,  the  FSCK  method  eliminates 
the  requirement  that  the  Planck  function,  or  solar  source  function,  be  constant 
over  the  spectral  intervals  to  be  integrated.  Consequently,  a  spectral  band  may 
be  arbitrarily  large  in  the  FSCK  approach,  even  encompassing  the  full  spectrum, 
as  long  as  the  assumption  of  correlation  in  the  gas  absorption  coefficients  between 
regions  of  the  atmosphere  remains  valid.  As  a  result  of  the  vertical  distribution 
of  absorbing  gases  in  the  atmosphere,  primarily  ozone  and  water  vapor,  the  cor¬ 
relation  assumption  breaks  down  when  the  full  shortwave  spectrum  is  treated  as 
a  single  band.  However,  the  breakdown  of  correlation  may  be  circumvented  by 
separating  the  spectrum  at  0.68  pm  into  two  bands,  one  band  in  which  ozone  con¬ 
tinuum  absorption  is  dominant  and  the  other  where  it  is  less  important  relative  to 
absorption  by  other  gases. 

Comparisons  of  two-band  FSCK  and  LBL  results  indicate  that  the  broadband 
FSCK  clear-sky  fluxes  and  heating  rates  are  accurate  to  better  than  1%  and  3%, 
respectively,  when  8192  quadrature  points  are  used  within  each  band.  Limiting 
the  first  band  from  0.24-0.68  pm  to  5  quadrature  points  and  the  second  band 
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from  0.68-4.60  /rm  to  10  quadrature  points,  the  accuracies  of  the  FSCK  fluxes  and 
heating  rates  become  approximately  1%  and  7%,  respectively.  For  low-  and  high- 
altitude  cloud  test  cases  we  found  that  gray  cloud  scattering  for  the  two  spectral 
bands  with  5-  and  10-point  quadratures  is  sufficient  to  reproduce  the  LBL  gener¬ 
ated  fluxes  and  heating  rates,  which  used  spectrally  varying  cloud  absorption  and 
scattering  properties,  to  approximately  1%  and  8%,  respectively.  Two  different 
approaches  to  modeling  absorption  by  cloud  drops  were  examined.  Explicitly  in¬ 
cluding  spectrally  varying  cloud  absorption  into  the  solar  source  function-weighted 
/^-distributions  resulted  in  realistic  in-cloud  heating  rates,  although  the  heating 
rates  were  consistently  less  than  LBL  calculated  rates.  Gray  cloud  absorption  pa¬ 
rameters  that  were  chosen  to  optimally  fit  LBL  results  for  one  atmosphere  but 
applied  to  different  atmospheres  or  cloud  combinations  also  closely  approximated 
LBL  in-cloud  heating  rates,  although  the  heating  rate  errors  were  less  consistent 
in  magnitude  and  sign  than  for  the  hrst  approach. 

The  FSCK  solution  requires  only  15  quadrature  points  per  layer,  a  40-90% 
reduction  in  computational  overhead  from  other  correlated  /^-distribution  models 
currently  in  use  (e.g.,  Fu  and  Liou  1992;  Mlawer  and  Clough  1997;  Mlawer 
and  Clough  1998;  Kato  et  al.  1999).  Hence,  the  FSCK  method  is  potentially 
important  to  numerical  modeling  of  radiative  transfer,  because  the  computational 
burden  of  the  radiative  transfer  calculations  is  strongly  linked  to  the  number  of 
spectral  bands  required.  The  next  step  in  the  application  of  the  FSCK  method  to 
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shortwave  atmospheric  radiative  transfer  is  to  develop  a  robust  radiative  transfer 
module  that  may  be  used  in  operational  numerical  weather  prediction  or  global 
climate  models.  Furthermore,  direct  comparisons  between  the  FSCK  method  and 
current  models  will  be  required  to  determine  the  relative  accuracy  of  the  models. 
By  explicitly  accounting  for  spectral  variability  of  the  solar  source  function,  the 
FSCK  method  has  the  potential  for  improving  accuracy  relative  to  traditional 
correlated  /^-distribution  models  which  assume  constant  solar  emission  over  each 
narrow  band. 


47 


Chapter  3 


Model  Development  for  General 
Atmospheres 

3.1  Design  Considerations 

To  advance  the  FSCK  model  discussed  in  Chapter  2  from  a  model  that  can  be  ap¬ 
plied  to  an  atmosphere  with  pre-specihed  pressures,  temperatures,  mixing  ratios, 
and  cloud  properties  to  an  operational  model  that  calculates  fluxes  and  heating 
rates  given  any  input  atmospheric  conditions,  there  are  many  factors  to  consider 
regarding  model  design.  For  the  study  discussed  in  Chapter  2  we  knew  the  prop¬ 
erties  of  each  atmospheric  layer  beforehand.  Therefore,  we  were  able  to  retrieve 
absorption  coefficients  for  mixtures  of  gases  with  pre-specihed  mixing  ratios  from 
LBLRTM  and  sort  them  to  produce  solar  source  function-weighted  cumulative  k- 
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distributions  for  each  model  layer.  Subsequently,  we  used  numerical  quadrature 
to  select  values  of  the  cumulative  /^-distributions  to  be  used  for  integration  of 
the  radiative  transfer  equation.  Since  production  of  absorption  coefficients  using 
LBLRTM,  sorting  them  into  /:-distributions,  and  summing  them  into  cumulative 
/^-distributions  are  computationally  expensive  procedures,  these  procedures  cannot 
be  performed  in  an  operational  model  at  runtime.  Therefore,  one  must  perform 
these  procedures  in  advance,  for  selected  pressures,  temperatures,  and  mixing  ra¬ 
tios,  storing  the  resulting  quadrature  values  of  absorption  coefficients  in  tables. 
When  the  operational  model  is  running,  interpolation  is  performed  between  stored 
values  to  absorption  coefficients  appropriate  for  the  actual  pressure,  temperature, 
and  mixing  ratio  values  in  each  atmospheric  layer. 

Important  considerations  for  implementing  these  procedures  are  the  size  and 
“shape”  of  the  stored  absorption  coefficient  table  and  the  method  of  combining 
the  optical  properties  of  constituent  species.  The  best  methods  for  interpolating 
stored  coefficients  to  the  actual  atmospheric  state  must  also  be  established  in  the 
operational  model.  Each  of  these  elements  must  be  examined  in  the  model  design 
stage  to  balance  accuracy  and  efficiency. 

An  example  of  the  importance  of  the  model  design  to  modern  Numerical 
Weather  Prediction  is  the  Weather  and  Research  Forecasting  (WRF)  model 
(Michalakes  et  al.  2001).  The  design  of  this  model  was  a  cooperative  venture 
amongst  a  number  of  major  national  and  educational  institutions.  As  such,  the 
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WRF  model  is  intended  to  be  the  next  generation  community  model  for  opera¬ 
tional  and  research  use.  With  such  a  broad  scope  of  applications,  strict  design 
standards  must  be  adhered  to  for  new  algorithms  implemented  into  the  model. 
One  such  standard  is  that  each  physical  parameterization  be  coded  within  a  single 
module.  In  the  case  of  our  full-spectrum  correlated  ^-distribution  radiative  trans¬ 
fer  algorithm,  the  implication  is  that  the  absorption  coefficient  table  cannot  be 
stored  in  an  external  data  table  to  be  read  by  the  parameterization  at  runtime  — 
the  table  must  be  coded  into  data  blocks  and  compiled  within  the  code.  The  size  of 
the  table,  which  is  now  important,  is  determined  by  the  number  of  pressures,  tem¬ 
peratures,  and  constituent  mixing  ratios  used  to  calculate  the  stored  absorption 
coefficients. 

The  “shape”  of  the  table  is  also  important,  where  “shape”  refers  to  the  specihc 
choice  of  pressures,  temperatures,  and  constituent  mixing  ratios  used  to  calculate 
the  stored  absorption  coefficients.  The  magnitude  of  the  separation  between  pres¬ 
sures,  temperatures,  and  constituent  mixing  ratios  contributes  to  the  quality  of  the 
interpolation  to  the  actual  atmospheric  state,  especially  for  absorption  coefficients 
not  linearly  related  to  the  atmospheric  state  variables.  Moreover,  the  range  of 
pressures,  temperatures,  and  constituent  mixing  ratios  must  encompass  as  many 
likely  natural  atmospheric  states  as  possible  to  avoid  extrapolation,  which  may 
introduce  large  errors.  As  a  result,  the  table  must  be  sufficiently  large  to  minimize 
heating  rate  errors,  yet  it  must  not  be  so  large  as  to  prohibit  its  inclusion  in  the 
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WRF  model. 


We  must  also  consider  the  method  of  combining  the  absorption  coefficients 
of  constituent  species.  The  LBLRTM  model  produces  absorption  coefficients  of 
individual  gases  in  spectral  space,  but  our  radiative  transfer  calculations  require 
absorption  coefficients  for  gas  mixtures  in  cumulative  ^-distribution,  or  g,  space. 
Adding  absorption  coefficients  of  constituent  species  in  spectral  space  is  simple 
and  accurate,  although  the  individual  gas  concentrations  are  not  known  in  ad¬ 
vance  within  an  operational  model  setting.  Because  orders  of  magnitude  fewer 
points  are  involved  in  g  space  compared  to  spectral  space,  it  may  be  possible  to 
combine  constituent  g'-space  values  within  the  run  time  of  the  numerical  model. 
Because  absorption  by  individual  species  may  be  correlated  and  one  must  account 
for  overlapping  absorption  lines  from  different  species,  combining  absorption  coef- 
hcients  in  ^'-space  requires  care  in  order  to  avoid  signihcant  errors.  Modest  and 
Riazzi  (2004)  address  this  topic  at  length,  and  we  will  discuss  it  further  in  Sections 
3.2  and  3.3. 

The  method  of  interpolation  must  also  be  considered  carefully,  balancing  accu¬ 
racy  against  efficiency.  Interpolating  from  a  table  of  stored  absorption  coefficients 
calculated  for  predetermined  pressures,  temperatures  and  constituent  mixing  ra¬ 
tios  requires  interpolation  in  at  least  three,  and  possibly  more,  dimensions.  While 
linear  interpolation  is  relatively  fast  compared  to  higher  order  interpolations  or 
spline  techniques,  the  relative  accuracies  of  the  different  methods  depends  on  the 
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problem  to  which  they  are  applied.  Finding  the  appropriate  method  for  a  particu¬ 
lar  problem  requires  many  “trial  and  error”  studies.  We  briefly  discuss  in  Sections 
3.2  and  3.3  the  motivation  for  the  interpolation  method  that  we  adopted  for  our 
FSCK  studies.  In  Section  3.4  we  describe  the  size  and  shape  of  the  correlated  k- 
coefficient  table  that  we  developed,  our  method  of  combining  constituent  radiative 
properties,  and  the  method  of  interpolation  that  we  chose  to  create  absorption  and 
scattering  coefficients  suitable  for  general  atmospheric  conditions. 

3.2  Precombined  Gases  with  Gray  Cloud 
Absorption 

In  Chapter  2  we  extracted  gas  constituent  properties  from  the  HITRAN  database 
and  combined  them  in  spectral  space  using  LBLRTM,  which  accounts  for  line 
overlap  and  continuum  absorption.  In  the  calculations  from  Chapter  2  we  cre¬ 
ated  constituent  radiative  properties  for  the  specihc  pressures,  temperatures,  and 
mixing  ratios  of  the  atmospheric  prohles  of  interest.  As  mentioned  in  Section  3.1, 
in  the  case  of  an  arbitrary  atmospheric  prohle  from  an  operational  model,  one 
does  not  know  the  values  of  the  atmospheric  state  variables  in  advance.  More¬ 
over,  the  computational  time  of  running  LBLRTM  and  combining  the  resulting 
gas  properties  is  far  too  expensive  during  the  runtime  of  an  operational  numerical 
model.  However,  we  can  tabulate  combined  constituent  A:-values  at  any  number 
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of  values  of  pressure,  temperature,  and  constituent  mixing  ratios  and  use  them  to 
interpolate  to  values  for  any  arbitrary  atmospheric  prohle. 

If  we  adopt  this  precombined  gas  property  approach,  we  must  compute  com¬ 
bined  ^-values  for  a  sufficient  number  of  combinations  of  pressures,  temperatures 
and  constituent  mixing  ratios  to  allow  for  reasonably  accurate  interpolated  absorp¬ 
tion  coefficients.  However,  too  many  tabulated  values  may  become  too  unwieldy 
to  use  in  practice.  For  example,  if  we  use  the  same  gases  that  we  used  in  Chapter  2 
and  store  15  double-precision  A;- values  for  each  combination  of,  say,  15  pressures,  15 
temperatures,  and  15  mixing  ratios  for  each  of  the  seven  gases,  the  resulting  table 
will  have  a  size  of  4.6  x  10^^  bytes.  Using  the  same  15  double-precision  quadrature 
points,  if  we  also  include  combinations  with  nongray  absorption  coefficients  for  15 
cloud-liquid  and  15  cloud-ice  drop-size  effective  radii,  the  table  increases  in  size  to 
1.0  X  10^®  bytes! 

One  way  to  reduce  the  size  of  the  table  is  to  treat  some  gases  as  minor  con¬ 
stituents,  and  to  set  their  mixing  ratio  to  a  single  “typical”  value  appropriate  for 
each  pressure  and  temperature  combination.  Major  constituents  are  those  gases 
with  multiple  mixing  ratios  for  each  pressure  and  temperature  combination.  Based 
on  the  preliminary  results  in  Chapter  2,  we  treat  water  vapor  as  a  major  constituent 
in  both  band  1  and  band  2,  ozone  as  a  major  constituent  in  band  1  and  a  minor 
constituent  in  band  2,  and  carbon  dioxide  as  a  minor  constituent  in  band  1  and 
a  major  constituent  in  band  2.  Minor  constituents  in  both  band  1  and  band  2 
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are  carbon  monoxide,  methane,  nitrous  oxide,  and  oxygen.  If  we  use  15  pressures, 
but  reduce  the  number  of  temperatures  to  10,  the  number  of  water  vapor  number 
concentrations  to  13,  the  number  of  ozone  mixing  ratios  in  band  1  to  6,  and  the 
number  of  carbon  dioxide  mixing  ratios  in  band  2  to  4,  while  further  limiting  the 
size  of  the  table  by  treating  cloud  absorption  as  gray,  i.e.,  making  a  separate  table 
for  cloud  optical  properties  that  is  independent  of  pressure  and  temperature  with 
one  value  per  band  for  each  lognormal  drop-size  distribution,  the  gas  absorption 
coefficient  table  is  reduced  to  a  more  reasonable  1.1  x  10®  bytes.  With  this  frame¬ 
work  in  mind  we  construct  Table  3.1  based  loosely  on  the  procedure  described 
by  Kato  et  ah  (1999).  Also  shown  in  Table  3.1  are  effective  radii  of  lognormal 
drop-size  distributions  of  liquid  and  ice  clouds  used  to  calculate  gray  cloud  scat¬ 
tering  properties  and  gray  or  nongray  cloud  absorption  independent  of  other  state 
variables. 

The  units  of  gas  concentration  used  to  tabulate  /^-values  for  different  gas  species 
requires  a  short  explanation.  Carbon  dioxide,  nitrous  oxide,  carbon  monoxide, 
methane,  and  oxygen  are  all  fairly  well-mixed  atmospheric  gases,  so  mixing  ratio 
(in  units  of  parts  per  million  by  volume,  [ppmv])  is  a  natural  unit  to  describe  their 
variation  through  the  atmospheric  prohle.  The  full  range  of  typical  mixing  ratio 
values  for  these  gases  can  be  resolved  adequately  with  only  a  few  values.  On  the 
other  hand,  water  vapor  concentrations  vary  signihcantly  through  the  atmospheric 
prohle,  so  many  values  of  mixing  ratio  are  necessary  to  adequately  resolve  the 
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full  range  of  water  vapor  mixing  ratios,  regardless  of  units  chosen.  We  would 
have  been  more  consistent  had  we  chosen  to  tabulate  water  vapor  mixing  ratios 
in  units  of  [ppmv].  However,  as  mentioned,  we  based  our  table  on  the  earlier 
work  of  Kato  et  ah  (1999),  in  which  water  vapor  is  tabulated  for  fixed  values  of 
number  concentration  (in  units  of  number  of  molecules  per  cubic  meter,  [m“^]). 
Since  number  concentration  is  the  most  appropriate  unit  for  radiative  transfer,  and 
in  fact  all  other  units  of  constituent  concentration  must  be  converted  to  number 
concentration  to  calculate  absorption  coefficients,  we  tabulate  water  vapor  A;- values 
for  fixed  values  of  number  concentration  in  units  of  [m“^],  thus  eliminating  the  need 
for  unit  conversion  at  model  runtime.  Ozone  is  neither  as  well  mixed  as  the  minor 
constituents,  nor  as  variable  as  water  vapor.  Consequently  we  tabulate  A;-values 
for  ozone  based  on  hxed  values  of  [ppmv],  but  use  more  mixing  ratios  than  the 
minor  constituents  to  resolve  the  typical  range  of  ozone  mixing  ratios. 

To  interpolate  from  stored  absorption  coefficients  associated  with  Table  3.1  to 
the  atmospheric  state  of  interest  requires  interpolation  in  four  dimensions,  i.e., 
pressure,  temperature,  water  vapor  number  concentration,  and  the  mixing  ratio  of 
either  ozone  for  band  1  or  carbon  dioxide  for  band  2. 

In  this  chapter  we  use  only  linear  interpolations  for  simplicity  and  efficiency, 
since  higher-order  interpolations,  while  likely  to  improve  accuracy,  are  certain  to 
reduce  efficiency,  especially  if  applied  to  this  four-dimensional  problem.  Similarly, 
cubic  spline  interpolations  would  likely  improve  accuracy,  but  cubic  spline  algo- 
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rithms  for  more  than  one  dimension  are  computationally  expensive,  as  well.  In  this 
chapter  we  determine  the  error  introduced  by  interpolating  tabulated  absorption 
coefficients  to  an  atmospheric  prohle  of  interest  by  applying  linear  interpolations. 
In  Chapter  4,  we  assess  whether  or  not  the  magnitude  of  the  interpolation  error  is 
acceptable  relative  to  state-of-the-art  operational  and  research  models.  If  it  is  not 
acceptable,  future  research  on  FSCK  atmospheric  applications  will  have  to  con¬ 
sider  implementing  a  higher-order  or  cubic-spline  interpolation  scheme,  or  perhaps 
a  hybrid  interpolation  scheme,  to  improve  accuracy. 

Note  that  a  four-dimensional  linear  interpolation  does  not  require  4  linear  in¬ 
terpolations  for  each  quadrature  point  of  each  layer;  rather,  it  requires  15  linear 
interpolations,  i.e.,  8  interpolations  in  pressure,  4  in  temperature,  2  in  water  vapor 
number  concentration,  and  1  in  the  mixing  ratio  of  the  other  major  gas. 

In  order  to  determine  the  error  introduced  by  interpolating  tabulated  absorp¬ 
tion  coefficients  to  an  atmospheric  prohle  of  interest  using  linear  interpolations, 
we  compare  the  heating  rates  resulting  from  the  use  of  our  table  of  precombined 
A:- values  to  the  heating  rates  produced  from  absorption  coefficients  computed  from 
the  precise  properties  of  the  atmospheric  prohle  of  interest,  which  we  call  the  “hard¬ 
wired”  atmosphere.  We  generate  2-band,  15  quadrature-point  solutions  for  both 
the  interpolated  and  “hardwired”  absorption  coefficients  and  use  the  “hardwired” 
results  as  the  standard  for  evaluating  the  results  based  on  interpolation.  This  ap¬ 
proach  allows  us  to  quantify  the  errors  produced  by  the  interpolation  procedure. 
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as  opposed  to  quantifying  the  loss  of  accuracy  relative  to  LBL  results  produced  by 
the  use  of  the  FSCK  algorithm  itself.  We  have  already  investigated  the  latter  type 
of  errors  in  Chapter  2. 

We  evaluate  our  interpolation  scheme  using  the  CLEAR,  CLOUD  A,  and 
CLOUD  B  cases  in  the  model  intercomparison  study  of  Barker  et  ah  (2003).  The 
atmospheric  thermodynamic  prohle  is  the  tropical  atmosphere  of  McClatchey  et  ah 
(1972)  with  62  layers  from  the  surface  to  approximately  70  km.  The  surface  albedo 
is  set  to  0.20.  In  the  CLOUD  A  case,  the  layer  between  3.5  km  and  4.0  km  contains 
a  nniform  overcast  clond  with  a  clond  water  mixing  ratio  of  Qc  =  0.159  g  kg“^.  The 
clond  water  is  modeled  as  lognormally  distribnted  pure  liquid  water  spheres  with 
an  effective  radius  of  refj  =  10  fim.  The  visible  optical  depth  that  results  is  ~10. 
The  CLOUD  B  case  uses  the  same  cloud-drop  size  distribution  as  for  the  CLOUD 
A  case,  but  hlls  the  layer  between  10.5  km  and  11.0  km  with  a  nniform  overcast 
clond  having  a  clond  water  mixing  ratio  of  Qc  =  0.034  g  kg“^  for  a  visible  optical 
depth  of 

As  the  results  in  Fig.  3.1a  demonstrate,  there  is  some  redaction  in  the  accuracy 
of  clear-sky  heating  rates  due  to  interpolation.  The  maximnm  clear-sky  heating 
rate  errors  relative  to  the  15-point  hardwired  solution  are  less  than  10%  when  the 
15-point  hardwired  solntion  heating  rate  has  a  value  greater  than  abont  1  K/day. 
Heating  rate  errors  don’t  exceed  5%  below  7  km.  Maximum  heating  rate  errors  are 
approximately  13%  at  the  top  of  the  atmospheric  prohle  (not  shown).  The  FSCK 


E 


D) 

0 

I 


0.0  0.5  1.0  1.5  0  2  4  6  8  10  0  2  4  6  8 

Heating  Rate  (K/day) 


Figure  3.1.  Broadband  heating  rates  (K/day)  for  the  FSCK  15-point  hardwired  solu¬ 
tion  (solid)  and  FSCK  calculations  using  precombined  gases  with  gray  cloud  absorption 
(dashed)  for  the  (a)  CLEAR,  (b)  CLOUD  A,  and  (c)  CLOUD  B  cases  of  Barker  et  ah 
(2003). 
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results  for  the  CLOUD  A  and  CLOUD  B  cases  are  less  promising  here,  as  compared 
to  the  results  in  Chapter  2,  when  gray  cloud  optical  properties  are  used  in  the 
simulations.  We  computed  gray  cloud  optical  properties  using  Mie  theory  applied 
to  the  Barker  et  ah  (2003)  cloud-drop  distributions.  Unfortunately,  we  could  not 
reduce  our  two-band  Mie  results  to  a  single  set  of  parameters  (i.e.,  extinction  cross 
section,  single-scattering  albedo,  and  asymmetry  parameter)  that  could  accurately 
reproduce  both  the  CLOUD  A  and  CLOUD  B  LBL  results.  We  could  tune  the 
gray  parameters  to  accurately  model  either  of  the  clouds,  but  parameters  tuned 
for  one  cloud  are  inappropriate  for  the  other.  The  implication  of  this  result  is 
that  any  operational  scheme  would  minimally  have  gray  cloud  parameters  that 
are  a  function  of  both  cloud-particle  size  and  their  vertical  location  within  the 
atmosphere.  Such  a  scheme  is  not  practical  to  implement.  Based  on  these  results 
we  do  not  develop  further  the  use  of  gray  cloud  optical  properties  for  the  two-band 
application  of  the  FSCK  method. 

3.3  Gases  Combined  On-the-Fly  with  Nongray 
Cloud  Absorption 

We  now  explore  the  advantages  and  disadvantages  of  another  approach  to  tabulat¬ 
ing  absorption  coefficients.  Rather  than  using  LBLRTM  to  combine  absorption  co¬ 
efficients  of  constituent  gases  in  spectral  space  prior  to  sorting  into  ^-distributions. 
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we  extract  and  sort  the  absorption  coefficients  of  each  individual  gas.  We  then  save 
the  /j-values  at  15  quadrature  points  for  each  gas  and  set  of  state  variables,  thereby 
creating  an  absorption  coefficient  table  for  each  gas.  We  also  create  nongray  ab¬ 
sorption  coefficient  tables  for  lognormally  distributed  cloud  liquid-  and  ice-particle 
size  distributions  with  a  range  of  effective  radii  in  the  same  manner  as  for  the  gases. 
We  subsequently  combine  the  gas  and  cloud  particle  absorption  coefficient  A:- values 
on-the-fly  using  the  method  of  Modest  and  Riazzi  (2004).  The  justihcation  for 
using  this  method  of  combining  nongray  cloud  absorption  and  gas  absorption  is  dis¬ 
cussed  in  Section  3.4.  Here,  we  use  two  different  implementations  of  this  method: 
we  combine  all  of  the  constituents  simultaneously  (Section  3.3.1)  and  we  combine 
the  constituents  serially  (Section  3.3.2). 

Because  the  constituents  properties  are  not  combined  for  various  mixing  ratios 
of  each  constituent,  the  tables  of  stored  coefficients  are  signihcantly  smaller  in 
size  for  this  approach.  The  A:-value  tables  for  each  constituent  are  at  most  three- 
dimensional,  with  axes  of  pressure,  temperature,  and  gas  mixing  ratio  given  in 
Table  3.1.  Because  the  tables  are  smaller,  we  treat  water  vapor,  ozone,  and  carbon 
dioxide  as  major  gases  in  both  band  1  and  band  2,  where  each  gas  is  indexed 
by  pressure,  temperature,  and  its  own  mixing  ratio.  We  treat  carbon  monoxide, 
methane,  nitrous  oxide,  and  oxygen  as  minor  constituents  for  both  bands,  so  their 
A:-value  tables  are  two-dimensional,  indexed  only  by  pressure  and  temperature. 
The  nongray  absorption  coefficient  tables  for  cloud-liquid  and  cloud-ice  particles 
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are  one-dimensional,  being  indexed  only  by  effective  radius.  We  use  these  same 
/j-value  tables  whether  we  combine  constituents  simultaneously  or  serially. 

As  an  alternate  implementation,  ozone,  carbon  dioxide,  and  the  other  minor 
gas  constituent  A:-values  could  be  calculated  for  a  single  pressure,  temperature,  and 
mixing  ratio  combination.  The  resulting  values  could  then  be  scaled  as  a  function 
of  pressure,  temperature,  and  mixing  ratio,  because  there  is  no  signihcant  self¬ 
broadening  associated  with  these  gases.  However,  although  scaling  of  A:-values  for 
these  gases  would  likely  be  more  efficient  than  interpolating  from  tabulated  values 
calculated  for  specihc  pressures,  temperatures,  and  mixing  ratios,  we  adhere  to  the 
current  implementation  to  avoid  introducing  additional  complexity  at  this  stage 
of  development.  We  leave  the  development  and  evaluation  of  scaling  functions  to 
subsequent  research. 

Because  the  /c-tables  for  the  individual  constituents  are  of  a  lower  dimension 
than  the  /c-table  for  precombined  gases,  interpolation  from  tabular  state  variables 
to  the  atmospheric  state  of  interest  is  also  performed  at  a  lower  dimension.  How¬ 
ever,  this  does  not  imply  fewer  interpolations.  In  the  case  of  precombined  gases 
with  gray  cloud  absorption,  we  performed  15  linear  interpolations  for  each  quadra¬ 
ture  point  at  each  atmospheric  layer.  Now  we  require  7  linear  interpolations  (i.e.,  4 
for  pressure,  2  for  temperature,  and  1  for  mixing  ratio)  for  each  major  constituent 
and  3  linear  interpolations  (i.e.,  two  for  pressure  and  1  for  temperature)  for  the 
remaining  minor  constituents.  Altogether,  we  now  have  a  total  of  24  linear  inter- 
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polations  for  each  quadrature  point  of  each  layer.  An  additional  interpolation  is 
necessary  for  layers  containing  either  a  liquid  or  an  ice  cloud,  while  2  additional 
interpolations  are  necessary  for  layers  with  a  mixed-phase  cloud. 

3.3.1  Constituents  Combined  Simultaneously 

The  method  that  we  chose  to  combine  tabulated  A:- values  obtained  from  constituent 
specie  cumulative  /^-distributions  was  adapted  from  Modest  and  Riazzi  (2004). 
The  derivation  of  the  Modest  and  Riazzi  (2004)  method  is  included  here  for  the 
convenience  of  the  reader. 

Taine  and  Souhani  (1999)  have  shown  that  the  product  of  the  transmissivities 
Ti  of  rii  individual  gas  species,  i.e.. 


rii 

rmix  =  n^*’  (3-1) 

i=l 

can  accurately  reproduce  direct  LBL  transmissivities  for  the  mixtnre.  This  is  true 
as  long  as  the  monochromatic  transmissivities,  and  hence  absorption  coefficients, 
of  the  constitnent  species  are  statistically  nncorrelated  across  the  spectral  interval 
for  which  the  Ti  are  computed.  Modest  and  Riazzi  (2004)  extended  Eq.  (3.1)  to 
a  method  of  mixing  cnmulative  /:-distributions  from  gas  species  with  nncorrelated 
absorption  coefficients.  Modest  and  Riazzi  (2004)  note  that  the  dehnition  of 
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transmissivity  in  terms  of  the  A:-distribntion  for  a  single  absorbing  specie, 


(3.2) 

is  also  the  dehnition  of  the  Laplace  transform  of  fi{k).  Using  Eqs.  (3.1)  and  (3.2), 
the  transmissivity  Tmix  of  a  mixture  of  /^-distributions  may  be  expressed  as  the 
product  of  the  Laplace  transforms  of  the  constituent  /^-distributions: 

Tbi 

Tmix  =  ^[/mix(/^)]  =  JJ (^)] •  (3.3) 

i=l 


In  terms  of  the  cumulative  /^-distribution  the  transmissivity  of  an  individual 
constituent  is 

e~’"^^dgi.  (3.4) 


For  a  mixture  of  two  constituents  this  becomes 


Tmix  =  C,[fram{k)]  =  6  ’^^^dgi  /  6 


'91=0  992=0 


where  ki  and  /c2  are  k-values  for  constituents  1  and  2,  respectively.  Modest  and 
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Riazzi  (2004)  use  the  integral  property  of  the  Laplace  transform  to  obtain 


C 


-1 


/„,,(«;)=  £[9„ixW1  =(  /  / 

'-^0  ^  ^Jgi=0Jg2=0 


-I  /■!  Q-{kl{9l)+k2{g2)]L 


(3.6) 


dgidgi- 


fgi=0  Jg2=0 


When  the  inverse  Laplace  transform  of  Eq.  (3.6)  is  taken  we  have 


fi^mix(^mix)  /  /  H[kniix  -  {h  +  k2)]dg2dgi,  (3.7) 

J gi=0  J 92=0 


where  /cmix  is  a  threshold  value  from  an  array  of  fc-bins  constructed  to  encompass 
the  full  range  of  /c-values  for  the  two  constituents  to  be  combined  in  a  model  layer 
and  dgi  and  dg2  are  the  quadrature  weights  associated  with  ki  and  k2-  In  Eq.  (3.7) 
H  is  the  Heaviside  step  function  given  by 

{1  for  ki  +  k2  <  fcmix 

(3.8) 

0  for  ki  +  k2>  /Cmix- 

The  function  g^ix  represents  the  cumulative  /o-distribution  of  the  mixture  to  be 
constructed  from  the  /c-values  for  each  of  the  constituent  gases. 

Modest  and  Riazzi  (2004)  construct  the  array  of  /c-bins  for  the  mixture,  /Cmix, 
using  a  power  distribution: 


A[(A;^ix)‘']  =  [(A;max)‘'  -  -  1),  (3.9) 
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and 


=  {krainy  +  *  A  [(/Cmix)'^]  (i  =  0,  •  ■  ■  ,  Nk^J,  (3.10) 

where  q  is  an  exponential  factor  set  to  0.095,  /Cmin  is  the  sum  of  the  minimum 
quadrature  values  of  the  two  constituents  in  the  model  layer,  and  k^ax  is  the  sum 
of  the  maximum  quadrature  values  of  the  two  constituents  in  the  model  layer.  The 
number  of  A:-bins,  is  important  to  the  accuracy  of  reconstructed  from 

the  values  of  ki  and  k2-  The  hner  the  resolution  in  /cmix,  the  more  precisely  the 
original  quadrature  values  will  be  captured  in  the  new  distribution. 

To  combine  i  constituents  Modest  and  Riazzi  (2004)  use  a  simple  generalization 
ofEq.  (3.7): 


^7mix(^mix) 


H\kxa\x  {ki  T  ■  ■  ■  T  ky^dgi  •  •  •  dgi. 


'31=0 


>gi=0 


(3.11) 


The  Heaviside  step  function  is  now  dehned  as 

{1  for  /Ci  H - h  /Cj  <  /Cmix 

(3.12) 

0  for  /ci  H - \-  ki>  kmix- 

We  still  compute  the  A;-bins  using  Eqs.  (3.9)  and  (3.10),  but  with  fcmin  and  k^ax 
obtained  from  summing  over  all  constituent  minimum  and  maximum  quadrature 
values,  respectively. 

As  a  result  of  the  nested  integrals  in  Eq.  (3.11),  the  innermost  loop  will  be 
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repeated  x  x  ■  ■  ■  x  times,  where  i  is  the  number  of  constituents 

and  Nk^  is  the  number  of  quadrature  points  for  the  constituent.  In  our  case, 
since  each  of  the  constituents  have  the  same  number  of  quadrature  points,  the 
innermost  loop  will  repeat  x  times  for  each  spectral  band,  where 

TVfcband  is  number  of  quadrature  points  for  each  spectral  band.  To  achieve 
reasonable  accuracy  we  choose  =  200.  For  5  quadrature  points  for  band 

1,  10  quadrature  points  for  band  2,  and  i  =  6  constituents,  which  includes  the 
three  major  gases,  the  minor  gases  treated  together  as  a  single  constituent,  and 
cloud  liquid-  and  ice-water  particles,  the  innermost  loop  will  repeat  approximately 
2.0  X  10®  times. 

In  our  application,  we  hrst  construct  a  new  cumulative  ^-distribution  for  the 
mixture,  fi'mix,  at  a  specihed  set  of  /cmix  values.  In  order  to  use  the  15  quadrature 
points  we  chose  in  Chapter  2  we  must  invert  the  array  5'mix(^mix)  to  /Cmix(fi'mix) 
and  interpolate  the  values  of  fc^ix  to  our  known  g'- values  (see  Table  2.1).  This 
procedure  requires  another  interpolation  for  each  quadrature  point  in  each  layer. 

The  heating  rates  for  the  Barker  et  ah  (2003)  atmospheric  prohles  resulting 
from  combining  gases  on-the-fly  simultaneously  with  nongray  cloud  absorption  are 
shown  in  Fig.  3.2.  The  heating  rate  errors  are  as  high  as  33%  for  the  CLEAR  case 
relative  to  the  15-point  hardwired  solution  in  the  lowest  10  km,  but  are  less  than 
10%  below  7  km.  This  approach  results  in  a  44%,  or  2.2  K/day,  in-cloud  error  in 
the  CLOUD  A  case  and  a  5%,  or  0.25  K/day,  in-cloud  error  in  the  CLOUD  B  case. 
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Figure  3.2.  Broadband  heating  rates  (K/day)  for  the  FSCK  15-point  hardwired  solution 
(solid)  and  FSCK  calculations  using  gases  and  nongray  cloud  absorption  combined  on- 
the-fly  simultaneously  (dashed)  for  the  (a)  CLEAR,  (b)  CLOUD  A,  and  (c)  CLOUD  B 
cases  of  Barker  et  al.  (2003). 


One  possible  source  of  error  present  in  this  approach,  which  is  not  present  in 
the  precombined  gas  approach,  is  in  the  interpolation  of  individual  gas  absorption 
coefficient  A:-values,  rather  than  the  mixture  absorption  coefficient  A:-values,  from 
tabulated  states  to  the  desired  atmospheric  state.  The  additional  interpolations 
necessary  in  the  process  of  combining  gases,  i.e.,  interpolating  values  of  kmw  from 
corresponding  gfmix-values  onto  our  15  quadrature  points,  may  introduce  additional 
errors.  Also,  it  is  possible  that  error  is  introduced  by  the  mixing  method  itself, 
because  it  assumes  absorption  coefficients  of  the  component  gases  are  statistically 
uncorrelated,  which  may  not  be  the  case.  Notwithstanding  all  of  these  errors,  in¬ 
cloud  absorption  is  more  accurate  with  the  current  treatment  of  combining  nongray 
cloud  absorption  coefficients  with  those  of  the  gases. 

3.3.2  Constituents  Combined  Serially 

A  slightly  different  implementation  of  the  approach  outlined  in  Section  3.3.1  is 
possible.  Rather  than  combining  the  constituents  simultaneously  through  Eq. 
(3.11),  we  simply  apply  Eq.  (3.7)  repeatedly.  That  is,  we  combine  the  hrst  two 
constituents  and  reassign  the  resulting  mixture  to  ki.  We  then  reapply  Eq.  (3.7) 
to  the  two-constituent  mixture  as  ki  and  the  ^-values  of  the  next  constituent  as 

k2. 

One  advantage  of  this  approach  is  computational  efficiency.  As  the  number 
of  constituents  increases,  the  computational  beneht  of  the  serial  approach  over 
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the  simultaneous  approach  also  increases.  Because  there  are  only  two  integrals 
in  Eq.  (3.7),  the  innermost  loop  will  be  repeated  (i  —  1)  x  x 

times,  where  again  i  is  the  number  of  constituents.  Because  of  the  improved 
computational  efficiency  of  this  method,  we  increase  =  2000  to  improve  the 

accuracy  of  the  combination  process.  For  our  six  constituents  (i.e.,  three  major 
gases,  the  minor  gases  treated  together  as  a  single  constituent,  and  cloud  liquid- 
and  ice- water  particles),  the  innermost  loop  will  repeat  only  1.2  x  10®  times  rather 
than  the  2.0  x  10®  times  required  for  simultaneous  combination  of  six  constituents. 

Another  advantage  of  this  approach  is  that  it  is  trivial  to  add  more  constituents. 
As  a  result,  adding  additional  gases  as  major  constituents  or  even  aerosols,  re¬ 
quires  only  the  addition  of  the  new  constituent  ^-tables  to  the  iteration  process. 
The  flexibility  of  this  approach  would  certainly  be  advantageous  in  the  design  of  a 
radiative  transfer  module  that  can  fully  interact  with  an  online  atmospheric  chem¬ 
istry  model.  A  disadvantage  is  an  increase  in  the  number  of  interpolations.  Each 
time  two  constituents  are  combined,  the  resulting  gmiY^miY)  must  be  inverted  and 
^mix(fl'mix)  must  be  interpolated  to  our  15  quadrature-point  values. 

The  heating  rates  resulting  from  combining  gases  on-the-fly  serially  with 
nongray  cloud  absorption  are  illustrated  in  Fig.  3.3.  The  heating  rate  errors 
are  less  than  26%  relative  to  the  15-point  hardwired  solution  for  the  CLEAR  case 
in  the  lowest  10  km,  with  a  34%  in-cloud  error  in  the  CLOUD  A  case  and  a  7% 
in-cloud  error  in  the  CLOUD  B  case.  Differences  in  the  results  between  the  se- 
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Figure  3.3.  Broadband  heating  rates  (K/day)  for  the  FSCK  15-point  hardwired  solution 
(solid)  and  FSCK  calculations  using  gases  and  nongray  cloud  absorption  combined  on- 
the-fly  serially  (dashed)  for  the  (a)  CLEAR,  (b)  CLOUD  A,  and  (c)  CLOUD  B  cases  of 
Barker  et  al.  (2003). 


rial  and  the  simultaneous  combination  approaches  are  numerical  in  origin.  The 
bin  resolutions  of  /cmix  are  different  and  the  threshold  values  of  /cmix  in  Eqs.  (3.7) 
and  (3.11)  are  met  in  a  different  sequence  in  the  two  approaches.  Another  reason 
for  differences  between  the  serial-  and  simultaneous-combination  approaches  is  the 
varying  number  of  interpolations  to  the  hnal  15  quadrature  points. 


3.4  The  Way  Forward 


Of  the  three  approaches  discussed  in  this  chapter,  the  clear-sky  heating  rate  errors 
were  the  smallest  for  the  precombined  gas  approach.  Whether  this  is  the  result  of 
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explicitly  accounting  for  line  overlap  by  combining  the  absorption  coefficients  in 
spectral  space  prior  to  sorting  or  the  fewer  interpolations  necessary  to  use  tabulated 
values  for  a  specihc  model  layer,  this  approach  appears  to  be  the  most  accurate 
for  clear-sky  calculations.  However,  when  we  used  this  method  with  gray  cloud 
properties,  in-cloud  heating  rate  errors  exceeded  125%.  As  mentioned  in  Section 
3.2,  if  we  include  nongray  cloud  absorption  in  the  precombined  mixtures,  the  size 
of  the  precombined  fc-table  will  become  unmanageable. 

A  compromise  approach  is  one  that  uses  precombined  gases,  but  combines 
nongray  cloud  absorption  on-the-fly.  In  this  approach  we  use  the  same  clear-sky 
table  developed  in  Section  3.2,  but  we  add  nongray  cloud  absorption  in  cloudy 
atmospheric  layers  serially  using  Eq.  (3.7). 

An  underlying  assumption  of  the  method  we  chose  to  combine  constituents, 
described  in  Section  3.3.1,  is  that  the  constituents  must  have  absorption  coeffi¬ 
cients  that  are  statistically  uncorrelated.  However,  the  absorption  coefficients  of 
cloud  drops  and  atmospheric  gases,  including  water  vapor,  may  not  be  entirely 
uncorrelated.  Therefore,  there  is  no  clear  theoretical  basis  for  applying  Eq.  (3.7) 
to  the  combination  of  cloud  absorption  coefficients.  To  determine  whether  or  not 
the  method  embodied  in  Eq.  (3.7)  is  appropriate  for  this  purpose  we  compare 
FSCK  band  2  transmissions  for  a  non-scattering  layer  of  a  clear-sky  gas  mixture, 
an  absorbing,  but  not  scattering,  cloud  with  a  lognormal  distribution  of  liquid 
cloud  drops  with  effective  radius  of  10  pm  and  a  varying  cloud-drop  number  con- 
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centration,  a  cloud-gas  mixture  whose  absorption  coefficients  have  been  combined 
monochromatically  and  the  product  of  the  clear-sky  gas  mixture  and  cloud  band 
2  transmissions  (Fig.  3.4b).  The  near  infrared  (wavelengths  1.0  /rm  and  3.0  /xm) 
cloud  absorption  optical  depths  resulting  from  the  variation  of  cloud-drop  number 
concentration  are  shown  in  Fig.  3.4a,  and  the  percent  error  between  the  product 
of  clear-sky  and  cloud  transmissions  relative  to  the  transmission  of  the  clear-  and 
clondy-sky  mixture  is  shown  in  Fig.  3.4c. 

As  Fig.  3.4c  illustrates,  errors  in  the  prodnct  of  constituent  transmissions  are 
less  than  10%  for  concentrations  less  than  86  cm“^  and  remain  less  than  15%  for 
concentrations  less  than  1000  cm“^.  We  expect  these  errors  to  be  representative 
of  other  cloud-drop  effective  radii,  as  well  as  for  ice  clouds  and  clouds  embedded 
in  air  with  slightly  different  gas  concentrations.  Based  on  the  magnitnde  of  these 
errors  we  proceed  to  use  the  method  of  combination  described  in  Section  3.3.1. 

The  heating  rates  resulting  from  the  use  of  precombined  gases  with  nongray 
clond  absorption  added  on-the-fly  serially  are  shown  in  Fig.  3.5.  The  resulting 
errors  are  less  than  13%  for  the  CLEAR  case  relative  to  the  15-point  hardwired 
solution  when  the  15-point  hardwired  solution  heating  rate  is  greater  than  about 
1  K/day.  When  the  15-point  hardwired  solntion  heating  rate  is  less  than  about 
1  K/day,  the  absolute  error  does  not  exceed  0.1  K/day.  The  in-clond  error  in  the 
CLOUD  A  case  is  25%  and  is  6%  in  the  CLOUD  B  case.  This  approach  is  the  one 
that  we  use  in  comparison  of  the  FSCK  method  with  state-of-the-art  operational 
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Figure  3.4.  (a)  Near  infrared  (wavelength  1.0  /xm  dotted,  and  wavelength  3.0  /xm  solid) 
absorption  optical  depths  for  a  lognormal  distribution  of  liquid-water  cloud  drops  with 
effective  radius  of  10  /um  as  a  function  of  number  concentration,  (b)  Broadband  near 
infrared  transmission  for  a  non-scattering  layer  consisting  of  a  clear-sky  gas  mixture 
(dotted),  a  lognormal  distribution  of  liquid- water  cloud  drops  with  effective  radius  of 
10  fim  (solid),  and  a  mixed  gas-cloud  layer  (dashed)  as  a  function  of  number  concentra¬ 
tion.  The  product  of  the  gas-mixture  and  cloud  broadband  transmissions  is  represented 
by  the  dash-dot  line,  (c)  Percent  error  of  the  product  of  the  clear-sky  gas  mixture  and 
cloud  broadband  transmissions  relative  to  the  broadband  transmission  for  the  mixed 
gas-cloud  layer  as  a  function  of  number  concentration. 
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Figure  3.5.  Broadband  heating  rates  (K/day)  for  the  FSCK  15-point  hardwired  solution 
(solid)  and  FSCK  calculations  using  precombined  gases  with  nongray  cloud  absorption 
combined  on-the-fly  serially  (dashed)  for  the  (a)  CLEAR,  (b)  CLOUD  A,  and  (c)  CLOUD 
B  cases  of  Barker  et  al.  (2003). 

and  research  radiative  transfer  methods  in  the  next  chapter. 
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Chapter  4 


Evaluation  versus  Operational 
and  Research  Models 


Based  on  the  results  from  Chapter  3  we  chose  to  implement  an  operational  FSCK 
approach  that  uses  a  table  of  quadrature  values  for  a  precombined  gas  mixture 
together  with  on-the-fly  serial  addition  of  nongray  cloud  absorption.  To  assess 
the  performance  of  the  two-band  FSCK  model  relative  to  other  currently  available 
operational  models,  we  now  compare  it  to  other  operational  GCM  or  NWP  model 
shortwave  radiative  transfer  parameterizations. 

A  number  of  shortwave  radiative  transfer  model  comparisons  have  been  under¬ 
taken  in  recent  years,  most  notably  the  Inter-comparison  of  Radiation  Codes  in 
Climate  Models  (ICRCCM)  program  of  the  World  Climate  Research  Programme 
(1984).  An  important  follow-up  study  was  performed  by  Barker  et  al.  (2003), 
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hereafter  referred  to  as  “the  ICRCCM-III  study.”  The  ICRCCM-III  study  in¬ 
cluded  results  from  25  1-D  models  and  compared  the  results  of  those  models  to 
calculations  of  four  3-D  Monte  Carlo  algorithms.  Since  the  3-D  models  were  not 
yet  able  to  perform  LBL  calculations,  the  four  3-D  models  were  also  compared 
to  results  calculated  by  the  CHARTS  model  (Clough  et  ah  1989),  a  1-D  LBL 
benchmark.  A  comparison  of  our  results  with  those  of  the  ICRCCM-III  study  will 
allow  us  to  characterize  the  performance  of  the  two-band  FSCK  algorithm  relative 
to  many  other  operational  and  research  shortwave  radiative  transfer  models  that 
are  currently  available.  This  comparison  will  be  discussed  in  more  detail  in  Section 
4.1. 

Although  many  of  the  models  included  in  the  ICRCCM-III  study  are  two- 
stream  solvers  similar  to  the  one  used  in  our  model  development  process 
(Zdunkowski  et  ah  1982),  differences  between  our  calculated  results  and  the 
ICRCCM-III  results  may  still  be  attributed  to  differences  in  specihc  model  solvers 
rather  than  the  treatments  of  gas  and  cloud  optical  properties,  which  is  the  focus  of 
this  work.  To  eliminate  differences  resulting  from  the  radiative  transfer  solver,  we 
also  compare  results  of  the  two-band  FSCK  approach  to  results  from  six  popular 
operational  and  research  GCM  and  NWP  model  radiative  transfer  algorithms  that 
we  have  embedded,  when  possible,  within  the  Zdunkowski  et  al.  (1982)  two-stream 
solver.  These  comparisons  allow  us  to  isolate  those  differences  that  result  from  dif¬ 
ferent  treatments  of  gas  and  cloud  optical  properties,  as  well  as  the  databases  upon 
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which  they  are  built.  We  describe  the  six  models  in  Section  4.2,  with  results  of 
our  comparisons  to  those  models  following  in  Section  4.3. 

4.1  Comparison  to  ICRCCM-III  Models 

In  onr  comparisons  of  the  performance  of  the  two-band  FSCK  model  to  the 
ICRCCM-III  operational  and  research  models,  we  adhere  as  closely  as  is  possi¬ 
ble  to  the  atmospheric  inpnts  for  each  case  study.  However,  there  are  technical 
featnres  of  onr  two-band  FSCK  approach  that  are  not  amenable  to  adaptation  to 
some  of  the  ICRCCM-III  analyses.  Therefore,  there  are  differences  in  implemen¬ 
tation  that  will  prodnce  differences  in  model  resnlts.  We  identify  differences  in 
specihc  resnlts  with  specihc  implementation  differences  whenever  possible. 

The  ICRCCM-III  stndy  compared  several  aspects  of  shortwave  radiative  trans¬ 
fer  parameterizations,  inclnding  treatments  of  clond  overlap  and  the  performance 
of  different  1-D  solvers  relative  to  the  3-D  solvers.  While  these  are  important 
considerations  in  the  development  of  a  complete  radiative  transfer  model,  in  this 
work  we  wish  to  focus  on  the  treatment  of  gas  and  clond  optical  properties.  There¬ 
fore,  we  perform  calculations  for  the  CLEAR  case  and  when  clonds  are  present, 
we  assnme  cloud  layers  are  simply  plane-parallel  and  homogeneous,  i.e.,  uniform 
overcast.  These  cases  constitute  a  subset  of  the  ICRCCM-III  study,  allowing  a  di¬ 
rect  comparison  to  the  ICRCCM-III  results  without  the  need  for  introdncing  any 
additional  complications. 
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Each  of  the  ICRCCM-III  cases  uses  the  lowest  62  layers  of  the  McClatchey 
et  ah  (1972)  tropical  atmosphere,  which  extends  from  the  surface  to  ~70  km.  A 
Lambertian  surface  with  gray  albedo  of  0.20  was  used.  The  pressure,  temperature, 
water  vapor  and  ozone  mixing  ratios  are  specihed  in  the  atmospheric  prohle,  with 
Barker  et  al.  (2003)  also  specifying  a  hxed  carbon  dioxide  mixing  ratio  of  360  ppmv. 
Barker  et  al.  (2003)  also  suggested  that  oxygen  be  included  in  the  ICRCCM- 
III  calculations,  but  they  did  not  provide  a  specihc  oxygen  mixing  ratio.  We 
include  in  our  model,  and  in  the  LBL  benchmark,  oxygen  mixing  ratios  adapted 
from  a  METEO  597A  class  project  at  the  Pennsylvania  State  University  (Fig. 
2.6).  Because  methane,  carbon  monoxide,  and  nitrous  oxide  are  built  into  the  gas 
mixtures  that  we  used  to  construct  the  precombined  gas  /s-table  for  our  operational 
FSCK  approach,  we  include  these  gases  as  minor  constituents.  The  inclusion  of 
these  gases  may  cause  some  difference  between  our  calculations  and  those  of  the 
ICRCCM-III  study. 

As  in  Chapter  3,  we  apply  the  operational  FSCK  approach  to  the  ICRCCM-III 
CLEAR,  CLOUD  A,  and  CLOUD  B  cases.  Cloud  properties  are  calculated  for 
lognormally  distributed  liquid- water  spheres  with  an  effective  radius  reff  of  10  /rm 
and  a  logarithmic  distribution  width  of  0.38.  In  the  CLOUD  A  case,  the  overcast 
cloud  hlls  the  layer  between  3.5  km  and  4.0  km  with  a  liquid  cloud  water  mixing 
ratio  of  Qc  =  0.159  g  kg“^,  resulting  in  a  visible  optical  depth  of  ~10  in  the  cloud 
layer.  In  the  CLOUD  B  case,  the  overcast  cloud  hlls  the  layer  between  10.5  km 
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and  11.0  km  with  =  0.034  g  kg“^,  resulting  in  an  in-cloud  visible  optical  depth 
of  ~1. 

In  the  comparison  below,  we  replicate  Figures  4,  5,  and  6  of  the  ICRCCM- 
III  study  (Barker  et  al.  2003)  using  model  outputs  both  from  the  Barker  et  ah 
(2003)  study  and  from  our  operational  two-band  FSCK  model  with  its  supporting 
LBL  calculations.  These  figures  display  results  of  the  1-D  models  as  a  function  of 
spectral  band,  as  Barker  et  al.  (2003)  have  divided  the  spectrum  into  two  bands, 
one  from  0.2-0. 7  /xm  and  the  other  from  0. 7-5.0  /xm.  As  described  in  Chapter  2, 
our  spectral  band  limits  are  0.24-0.68  /rm  and  the  other  from  0.68-4.60  /rm.  The 
full  spectrum  results  are  also  shown,  where  the  limits  of  the  ICRCCM-III  models 
are  0. 2-5.0  /xm  and  the  limits  of  our  model  are  0.24-4.60  /rm.  The  difference 
between  the  spectral  limits  is  one  source  of  differences  in  the  spectral  band  results. 

The  original  ICRCCM-III  hgures  also  show  results  for  one  of  the  3-D  Monte 
Carlo  codes  and  CHARTS  LBL  benchmark  calculations.  Since,  for  this  study, 
we  are  not  interested  in  differences  resulting  from  different  solvers,  we  do  not 
include  the  Monte  Carlo  and  CHARTS  results  in  our  hgures.  For  the  benchmark 
calculations  we  use  instead  our  own  LBL  calculations  with  the  Zdunkowski  et  al. 
(1982)  two-stream  solver.  Finally,  while  we  have  included  results  of  the  1-D  models 
as  a  function  of  cosine  of  solar  zenith  angle,  /xq,  we  generated  the  results  for  the 
FSCK  model  only  for  /io  =  0.5. 

Figure  4.1a-c  shows  top-of-atmosphere  (TOA)  albedo,  surface  absorptance,  and 
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Figure  4.1.  (a)  Top-of-atmosphere  albedo,  (b)  surface  absorptance,  and  (c)  atmo¬ 

spheric  absorptance  for  the  CLEAR  case  as  a  function  of  cosine  of  solar  zenith  angle 
ICRCCM-III  1-D  model  median  values  are  represented  by  solid  lines,  while  the  upper 
and  lower  quartiles  are  represented  by  dashed  lines.  Our  LBL  calculations  are  plotted  as 
open  circles  and  FSCK  calculations  are  plotted  as  solid  circles  for  ^o=0.5.  (d)  Heating 
rate  profiles  for  the  CLEAR  case.  Here,  the  thick  solid  line  represents  our  LBL  bench¬ 
mark  calculations,  whereas  the  thin  solid  and  dashed  lines  are  the  same  as  for  (a)-(c). 
The  FSCK  calculations  are  plotted  as  a  dotted  line. 
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atmospheric  absorptance  for  the  CLEAR  case  as  a  function  of  po-  In  fhe  figures 
the  median  values  of  the  1-D  ICRCCM-III  models  along  with  the  25%  and  75% 
percentiles  are  illustrated,  as  well  as  the  FSCK  calculations  for  po  =  0.5.  The 
FSCK  broadband  TOA  albedo  agrees  to  within  0.001  of  the  median  1-D  value, 
but  the  FSCK  atmospheric  absorptance  is  high  and  the  surface  absorptance  is 
low  relative  to  the  median  values  of  the  ICRCCM-III  1-D  models;  this  difference 
is  also  apparent  in  the  near  IR  band  2  results.  FSCK  calculations  for  both  the 
broadband  and  band  2  surface  and  atmospheric  absorptances  are  closer  to  the 
LBL  calculations  than  are  the  1-D  median  values.  Barker  et  ah  (2003)  noted  that 
the  1-D  codes  in  the  ICRCCM-III  study  markedly  underestimated  atmospheric 
absorptance  and  overestimated  surface  absorptance  relative  to  the  Monte  Carlo 
and  LBL  benchmarks,  and  postulated  that  the  underestimations  were  likely  due 
to  lack  of  water  vapor  continuum  absorption  and  parameterizations  based  on  old 
spectroscopic  databases  relative  to  the  newer  parameterizations  used  by  the  3-D 
and  benchmark  codes.  Since  the  FSCK  model  includes  water  vapor  continuum 
absorption  and  is  based  on  the  same  fairly  recent  spectroscopic  data  as  our  LBL 
benchmark,  i.e.,  the  HITRAN-96-JPL  database  (Rothman  et  ah  1998),  we  ex¬ 
pect  the  FSCK  model  to  have  somewhat  higher  atmospheric  absorptance  than  the 
median  of  the  ICRCCM-III  1-D  models. 

The  heating  rates  for  the  CLEAR  case  are  shown  in  Fig.  4.  Id.  Although  we 
do  not  use  the  same  LBL  benchmark  as  employed  in  the  ICRCCM-III  case,  our 
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benchmark  calculations  demonstrate  the  same  feature  remarked  upon  by  Barker  et 
al.  (2003).  Below  ~  350  hPa  “the  benchmark  values  are  systematically  larger  than 
the  vast  majority  of  ID  estimates.”  Barker  et  al.  (2003)  again  attributed  those 
differences  to  poor,  or  old,  water  vapor  line  and  continuum  absorption  information 
used  by  the  1-D  models.  In  the  case  of  the  FSCK  model,  in  which  we  use  the 
same  water  vapor  line  and  continuum  absorption  as  our  LBL  benchmark,  our 
heating  rates  are  generally  closer  to  the  benchmark  calculations  than  the  other 
1-D  models.  The  mean  bias  heating  rate  error  relative  to  the  LBL  benchmark  over 
all  atmospheric  layers  is  -0.27  K/day  for  the  median  of  the  1-D  models,  but  is  only 
-1-0.04  K/day  for  the  FSCK  model.  Another  likely  source  for  the  improvement 
of  FSCK  over  the  ICRCCM-III  1-D  models  is  the  fact  that  the  FSCK  approach 
explicitly  accounts  for  the  spectral  variation  of  the  solar  irradiance.  Where  the 
FSCK  heating  rates  differ  from  the  LBL  calculations,  errors  are  primarily  due 
to  interpolations  from  tabulated  absorption  coefficients  as  discussed  in  Chapter 
3,  rather  than  a  dehciency  in  the  FSCK  algorithm  itself.  Recall  that  the  FSCK 
algorithm  produced  more  accurate  results  relative  to  LBL  calculations  when  the 
exact  atmospheric  prohle  was  prespecihed  as  in  Chapter  2. 

Figures  4.2  and  4.3  show  results  for  the  CLOUD  A  and  CLOUD  B  cases, 
respectively.  Note  that  we  have  reduced  the  vertical  extent  of  the  y-axes  in  Figs. 
4. 2d  and  4.3d  in  order  to  show  heating  rates  within  and  immediately  above  and 
below  the  cloud  layers.  As  in  Fig.  4.1,  Figs.  4.2a-c  and  4.3a-c  show  median  and 
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Figure  4.2.  (a)  Top-of-atmosphere  albedo,  (b)  surface  absorptance,  and  (c)  atmospheric 
absorptance  for  the  CLOUD  A  case  as  a  function  of  cosine  of  solar  zenith  angle  ^o- 
ICRCCM-III  1-D  median  values  are  represented  by  solid  lines,  while  the  upper  and 
lower  quartiles  are  represented  by  dashed  lines.  Our  LBL  calculations  are  plotted  as 
open  circles  and  FSCK  calculations  are  plotted  as  solid  circles  for  ^o=0.5.  (d)  Heating 
rate  profiles  for  the  CLOUD  A  case.  Here,  the  thick  solid  line  represents  our  LBL 
benchmark  calculations,  whereas  the  thin  solid  and  dashed  lines  are  the  same  as  for 
(a)-(c).  The  FSCK  calculations  are  plotted  as  a  dotted  line. 
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Figure  4.3.  (a)  Top-of-atmosphere  albedo,  (b)  surface  absorptance,  and  (c)  atmospheric 
absorptance  for  the  CLOUD  B  case  as  a  function  of  cosine  of  solar  zenith  angle  ^o- 
ICRCCM-III  1-D  median  values  are  represented  by  solid  lines,  while  the  upper  and 
lower  quartiles  are  represented  by  dashed  lines.  Our  LBL  calculations  are  plotted  as 
open  circles  and  FSCK  calculations  are  plotted  as  solid  circles  for  ^o=0.5.  (d)  Heating 
rate  profiles  for  the  CLOUD  B  case.  Here,  the  thick  solid  line  represents  our  LBL 
benchmark  calculations,  whereas  the  thin  solid  and  dashed  lines  are  the  same  as  for 
(a)-(c).  The  FSCK  calculations  are  plotted  as  a  dotted  line. 
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quartile  values  for  the  ICRCCM-III  1-D  models  rather  than  mean  values.  This 
explains  the  apparent  contradiction  of  FSCK  atmospheric  absorptance  values  in 
Fig.  4.2c  that  are  low  relative  to  the  1-D  model  median  values  for  both  band  1 
and  band  2,  but  are  virtually  identical  to  the  full  shortwave  1-D  median  value. 

In  both  cloud  cases  the  heating  rate  profiles  demonstrate  that  the  FSCK  treat¬ 
ment  of  cloud  absorption  as  nongray  is  more  accurate  relative  to  LBL  calculations 
than  the  ICRCCM-III  1-D  models,  which  tend  to  overabsorb  in  the  cloud  layers. 
In-cloud  heating  rate  errors  for  the  CLOUD  A  case  were  6.1%  for  the  FSCK  ap¬ 
proach  and  16.0%  for  the  1-D  model  median.  For  the  CLOUD  B  case  in-cloud 
heating  rate  errors  were  -3.5%  for  the  FSCK  approach  and  15.0%  for  the  1-D 
model  median.  The  most  likely  source  of  error  for  the  ICRCCM-III  1-D  models  is 
averaging  cloud  properties  over  spectral  intervals  where  the  properties  are  chang¬ 
ing  rapidly.  To  demonstrate  that  the  overestimate  of  heating  by  the  ICRCCM-III 
1-D  models  relative  to  our  LBL  calculations  is  not  the  result  of  the  ICRCCM-III 
models  extending  their  calculations  in  the  near  infrared  to  5.0  pm  versus  the  FSCK 
upper  limit  of  4.6  pm,  we  performed  LBL  calculations  from  a  wavelength  of  4.6- 
5.0  pm  and  added  the  absorbed  energy  to  our  LBL  results.  The  clear  sky  heating 
rates  increased  by  0.01  K/day  or  less,  while  the  in-cloud  heating  rates  increased 
by  approximately  0.09  K/day,  or  1.6%  for  CLOUD  A  and  by  0.13  K/day  or  2.5% 
for  CLOUD  B. 


Given  the  higher  in-cloud  heating  rates  of  the  ICRCCM-III  models  relative  to 


FSCK  and  LBL  calculations,  one  might  expect  ICRCCM-III  1-D  model  column 
atmospheric  absorptance  to  exceed  that  of  the  FSCK  calculations.  This  is  the  case 
for  bands  1  and  2  of  CLOUD  A  and  band  1  of  CLOUD  B.  However,  FSCK  atmo¬ 
spheric  absorptance  for  CLOUD  B  exceeds  that  of  the  ICRCCM-III  1-D  models 
for  band  2  and  for  the  full  shortwave.  We  attribute  this  result  to  poor  water  vapor 
line  and  continuum  absorption  in  the  ICRCCM-III  1-D  models,  just  as  we  did  for 
the  source  of  differences  in  the  CLEAR  case  comparisons.  In  both  of  the  cloud 
cases  for  all  of  the  spectral  bands  the  median  1-D  atmospheric  absorptance  is  less 
than  the  LBL  benchmark.  As  pointed  out  by  Barker  et  ah  (2003),  the  ICRCCM- 
III  1-D  models  “underestimate  systematically  relative  to... CHARTS  by  amounts 
similar  to  those  for”  the  CLEAR  case. 

As  one  would  expect,  with  less  cloud  absorption  the  FSCK  results  for  the 
CLOUD  A  and  CLOUD  B  cases  show  more  energy  absorbed  at  the  surface  relative 
to  the  ICRCCM-III  1-D  models.  Lower  FSCK  TOA  albedo  values  relative  to 
the  ICRCCM-III  results,  especially  for  the  CLOUD  A  case,  have  a  number  of 
possible  sources.  Primary  among  them  is  less  than  optimal  tuning  of  the  gray 
cloud  asymmetry  parameter  (Section  2.2.3).  This  might  also  contribute  to  FSCK 
atmospheric  absorptance  values  that  underestimate  those  calculated  by  LBL  for 
the  CLOUD  A  and  B  cases. 
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4.2  State-of-the-Art  Operational  and  Research 


Models 

While  the  comparisons  with  the  ICRCCM-III  cases  demonstrate  that  the  FSCK 
approach  produces  reasonable  results  relative  to  other  methods,  we  cannot  easily 
use  the  ICRCCM-III  study  to  separate  differences  between  the  solvers,  the  treat¬ 
ments  of  gas  and  cloud  optical  properties,  and  any  manipulations  of  the  input  fields. 
To  eliminate  the  solver  and  manipulations  of  the  input  fields  as  the  sources  of  any 
differences,  we  picked  six  NWP  and  research-grade  schemes  commonly  used  for 
handling  gas  and  cloud  absorption  and  scattering  and  we  embedded  them  within 
the  Zdunkowski  et  ah  (1982)  two-stream  solver.  We  can  draw  much  stronger 
conclusions  regarding  treatment  of  gas  absorption  using  the  results  from  these  six 
models,  as  well  as  our  FSCK  and  LBL  results,  all  of  which  use  the  same  solver. 

The  six  shortwave  radiative  transfer  algorithms  that  we  chose  for  this  study  are 
the  Dudhia  shortwave  scheme  (Dudhia  1989),  the  Geophysical  Fluid  Dynamics  Lab 
(GFDL)-Eta  scheme  (Lacis  and  Hansen  1974),  the  Community  Atmosphere  Model 
version  2  (CAM2)  and  version  3  (CAMS)  schemes  (Briegleb  1992)  counted  as  one 
model,  the  correlated  ^-distribution  model  developed  by  Fu  and  Liou  (1992),  the 
correlated  /^-distribution  model  developed  by  Kato  et  al.  (1999),  and  the  shortwave 
correlated  /^-distribution  model  Rapid  Radiative  Transfer  Model  (RRTM  SW; 
Mlawer  and  Clough  1997;  Mlawer  and  Clough  1998).  A  breakdown  for  these 


models  of  the  treatment  of  gas  and  cloud  optical  properties  by  spectral  interval 
is  depicted  in  Fig.  4.4.  In  the  hgure,  spectral  intervals  not  treated  by  a  model 
are  indicated  by  a  diagonal  line  through  the  interval.  The  hrst  three  models  rely 
heavily  on  the  empirical  methods  of  Lacis  and  Hansen  (1974),  while  the  other 
models,  separated  from  the  first  three  by  a  blank  line  in  Fig.  4.4,  are  three  different 
implementations  of  the  correlated  /^-distribution  approach. 

The  codes  for  the  models  used  in  this  comparison  were  obtained  from  a  variety 
of  sources,  including  from  modules  within  operational  NWP  models.  Because  our 
focus  here  is  on  optical  property  generation,  we  extracted  only  those  parts  of  the 
codes  necessary  for  our  comparison.  In  the  case  of  the  GFDL-Eta  scheme,  the  code 
for  cloudy  conditions  was  heavily  interwoven  through  a  number  of  subroutines  with 
dependencies  upon  other  modules.  We  felt  that  the  time  spent  extracting  the  cloud 
optical  property  code  for  the  GFDL-Eta  scheme  was  better  spent  elsewhere.  As  a 
consequence,  no  results  are  included  for  the  GFDL-Eta  scheme  for  the  GLOUD  A 
and  GLOUD  B  cases. 

The  work  by  Lacis  and  Hansen  (1974)  is  a  fundamental  source  for  three  of  the 
models  used  in  this  comparison,  although  each  of  the  three  models  uses  different 
aspects  of  the  work.  Because  of  the  relevance  of  their  work  to  current  operational 
models,  we  will  take  a  close  look  at  the  development  of  their  treatments  of  gas 
absorption  and  how  those  treatments  are  applied  in  the  three  models.  We  will  also 
show  how  a  number  of  the  six  models  in  our  study,  including  that  of  Lacis  and 
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Figure  4.4.  Schematic  diagram  of  spectral  bands  for  six  popular  operational  and  research  shortwave  radiative  transfer  models 
and  the  FSCK  model.  Diagonal  lines  indicate  spectral  regions  not  included  by  the  model. 


Hansen  (1974),  discuss  or  implement  ideas  that  are  similar  to  the  FSCK  method 
for  the  specihc  case  of  solar  radiation  and  point  out  where  those  models  differ  from 
the  FSCK  approach. 

4.2.1  Dudhia  Shortwave  Scheme 

The  Dudhia  shortwave  scheme  (Dudhia  1989),  which  is  used  in  the  hfth  gener¬ 
ation  Penn  State  University /National  Center  for  Atmospheric  Research  (NCAR) 
Mesoscale  Model  (MM5)  (Grell  et  ah  1994)  and  the  Weather  Research  and  Fore¬ 
casting  (WRF)  model  (Michalakes  et  ah  2001),  is  based  on  the  absorption  by 
water  vapor  in  the  near  infrared  according  to  Lacis  and  Hansen  (1974).  The  treat¬ 
ment  of  water  vapor  absorption  by  Lacis  and  Hansen  (1974)  is,  in  turn,  based  on 
work  by  Howard  et  al.  (1956)  and  Yamamoto  (1962).  Yamamoto  (1962)  dehned 
absorptivity  as  the  ratio  of  the  energy  absorbed  by  the  entire  vertical  air  column 
to  the  solar  constant,  assuming  normal  incidence.  Applying  that  dehnition,  Ya¬ 
mamoto  (1962)  took  the  measured  absorption  data  of  Howard  et  al.  (1956),  as 
well  as  data  for  two  additional  weak  absorption  bands  not  included  in  the  earlier 
measurements,  and  weighted  them  by  measured  incoming  solar  energy  in  the  same 
spectral  intervals  to  get  fractional  absorptivities.  Yamamoto  (1962)  then  summed 
the  fractional  absorptivities  to  get  a  total  absorptivity  for  the  near  infrared  (0.7- 
4.0  fim)  and  plotted  a  curve  of  total  absorptivity  as  a  function  of  precipitable  water 
vapor.  Lacis  and  Hansen  (1974)  then  £t  a  curve  using  a  least  squares  fit  method 
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to  the  curve  of  Yamamoto  (1962),  arriving  at  a  function  for  absorption  due  to 


water  vapor  given  by 


(1  +  141.5|/)0-635  +  5.925|/’ 

where  A^viv)  is  the  fraction  of  incident  solar  flux  that  is  absorbed  by  water  vapor, 
i.e.,  the  absorptivity,  as  a  function  of  the  water  vapor  amount,  y  (in  units  of 
centimeters  of  precipitable  water  vapor). 

Molecular  and  aerosol  scattering  in  the  Dudhia  model  are  treated  as  gray  and 
proportional  to  the  atmospheric  mass  path  length,  which  accounts  for  the  solar 
zenith  angle,  with  a  constant  giving  10%  reflectivity  for  1  atm.  Absorptivity  and 
reflectivity  due  to  clouds  are  provided  by  tabulated  functions  of  the  cosine  of 
the  zenith  angle,  po,  and  \n{w/fio),  where  w  [gm~‘^]  is  the  vertically  integrated 
cloud  water  path.  These  tabulated  functions  are  derived  from  theoretical  results  of 
Stephens  (1978).  The  gas  and  cloud  optical  properties  used  in  the  Dudhia  model 
are  not  easily  put  into  a  form,  i.e.,  gas  absorption  optical  depth,  gas  scattering 
optical  depth,  cloud  extinction  optical  depth,  cloud  asymmetry  parameter,  and 
cloud  single-scattering  albedo,  needed  by  the  Zdunkowski  et  al.  (1982)  two- 
stream  solver.  Consequently,  we  calculated  heating  rates  using  the  native  solver 
of  the  Dudhia  scheme  for  the  comparisons  in  Section  4.3. 

The  Dudhia  solver  simply  attenuates  the  downward  flux  layer  by  layer.  First, 
the  top  layer  is  treated  as  though  it  is  the  entire  atmospheric  column,  and  the  ab- 
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sorption  for  that  column  is  calculated.  Next,  the  atmospheric  thickness  is  extended 
to  include  the  top  two  layers,  and  the  absorption  for  that  thickness  is  calculated, 
and  so  on.  The  heating  rates  for  a  layer  are  determined  by  subtracting  absorption 
by  water  vapor  and  cloud  for  the  previously  calculated  thickness  from  that  of  the 
current  thickness.  The  simplicity  of  this  model  gives  it  a  signihcant  advantage 
over  other  models  included  in  this  comparison  with  respect  to  efficiency,  but  this 
efficiency  comes  at  the  expense  of  accuracy,  as  will  be  shown  in  Section  4.3. 

4.2.2  GFDL-Eta  Scheme 

The  GFDL-Eta  scheme  used  in  the  National  Centers  for  Environmental  Prediction 
(NCEP)  Eta  model  and  the  WRF  model  is  also  based  on  the  work  of  Lads  and 
Hansen  (1974),  adopting  from  that  work  the  parameterization  of  absorption  by 
ozone  and  by  water  vapor  in  cloudy  skies.  This  model  uses  two  spectral  bands.  One 
band,  covering  ultraviolet  (A  <  0.35  /rm)  and  visible  wavelengths  (0.5  <  0.7  /rm), 
accounts  for  absorption  by  ozone  and  water  vapor,  as  well  as  molecular  “Rayleigh” 
scattering.  A  second  band,  for  near  infrared  wavelengths,  accounts  for  absorption 
by  water  vapor  and  carbon  dioxide. 

Using  a  least  square  ht  to  solar  flux-weighted  measured  ozone  absorptivities, 
Lacis  and  Hansen  (1974)  express  the  fraction  of  solar  flux  across  the  ultraviolet 
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wavelengths  absorbed  by  the  Hartley  and  Huggins  ozone  bands  as 


1.082a; 


X  = 


+ 


0.0658a; 


(1  +  138.6a;)0-805  (1  +  103.6a;)3  ’ 


(4.2) 


where  x  is  ozone  amount  (in  units  of  cm  at  293.16  K  and  1013.25  hPa).  The 
fraction  of  solar  flux  across  visible  wavelengths  absorbed  by  the  Chappuis  ozone 
band  is  given  by 


AT/(x) 


0.02118a; 

(1  +  0.042a;  +  0.000323a;2)  ■ 


(4.3) 


Total  ozone  absorption  then  is  given  by 


aux)  =  a::(x)+a:-(x). 


(4.4) 


Absorption  by  carbon  dioxide  in  the  near  infrared  is  parameterized  according 
to  Sasamori  et  al.  (1972),  who  use  the  relation 

a  =  2.35  X  10-=^(m  +  O.Om)^-^^^  -  7.5  x  10-^  (4.5) 


where  a  is  the  fractional  absorptivity  and  u  is  the  pressure-corrected  path  length 
of  carbon  dioxide  (in  units  of  cm  reduced  to  STP).  This  relation  is  an  empirical 
curve  £t  to  the  data  of  Burch  et  al.  (1960). 

The  GFDL-Eta  model  uses,  for  clear  and  cloudy  conditions,  an  approach  for 
absorption  by  water  vapor  in  the  near  infrared  that  is  adapted  from  the  cloudy-sky 


94 


formulation  by  Lacis  and  Hansen  (1974),  who  use  a  “^-distribution”  consisting 
of  8  points  derived  from  measured  absorptivities.  Lacis  and  Hansen  (1974)  did 
not  have  a  “database  of  monochromatic  absorption  coefficients  and  their  pressure- 
temperature  dependence.”  So,  rather  than  construct  a  true  /^-distribution,  they 
used  a  least-square  search  to  numerically  invert  the  equation 

poo  ^ 

A^^{y)  =  l-  p{k)e~^^dk  ^  1 -^p{kn)e~^”^Akn:  (4.6) 

Jo 


to  obtain  N  values  for  the  discrete  probability  distribution  p{kn)  and  the  corre¬ 
sponding  absorption  coefficient  kn-  The  absorption  for  water  vapor  is  corrected  for 
pressure  and  temperature  by  using  an  effective  absorber  amount  given  by 


Ves  =  y 


(4.7) 


This  approach  is  effectively  the  weighted-sum-of-gray-gases  method  (Hottel  and 
Sarohm  1967;  Modest  1991).  Lacis  and  Hansen  (1974)  state  that  “the  expression 
p{k)dk  is  the  fraction  of  the  incident  flux  that  is  associated  with  an  absorption 
coefficient  between  k  and  k  -|-  dk”  in  essence  describing  the  FSCK  method  for  the 
special  case  of  shortwave  radiation.  Had  a  “database  of  monochromatic  absorption 
coefficients  and  their  pressure-temperature  dependence”  been  available  to  them, 
Lacis  and  Hansen  (1974)  might  have  developed  the  FSCK  method  of  Modest 
(2003). 
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The  weights,  p{kn),  used  in  the  GFDL-Eta  model  have  been  altered  from  those 
given  by  Lacis  and  Hansen  (1974).  Moreover,  three  additional  A:- values  have 
been  added  to  improve  the  accuracy  of  absorption  by  water  vapor  (Y.-T.  Hon 
2004,  personal  communication).  Rayleigh  scattering  is  accounted  for  as  Lacis 
and  Hansen  (1974)  describe  for  clear-sky  conditions  and  is  only  calculated  for  the 
ultraviolet-visible  band  and  only  for  the  lower  atmosphere.  The  reflection  function 
for  the  lower  atmosphere  for  the  UV-VIS  band  is  given  by 


R{fj,o) 


0.219 

1  -|-  0.816/io  ’ 


(4.8) 


where  /io  is  the  cosine  of  the  solar  zenith  angle.  This  formula  is  again  obtained 
by  an  iterative  least-square  ht  to  obtain  the  best  agreement  between  the  vertical 
absorption  prohles  for  their  parametric  formula  and  their  benchmark. 

As  with  the  Duhdia  scheme,  the  gas  absorption  and  scattering  properties  are 
hrmly  integrated  with  the  solver,  an  adding  and  doubling  scheme  also  described 
by  Lacis  and  Hansen  (1974).  So  rather  than  attempt  to  disentangle  the  properties 
from  the  solver  and  implement  them  into  the  Zdunkowski  et  al.  (1982)  two-stream 
solver,  we  use  heating  rates  calculated  by  this  model’s  solver  for  our  comparisons 
in  Section  4.3. 


96 


4.2.3  CAM2  and  CAMS  Schemes 


The  Community  Atmosphere  Model  version  2  (CAM2;  Collins  1998),  which  is  the 
NCAR  community  climate  model,  includes  a  radiative  transfer  parameterization 
that  is  an  updated  version  of  the  parameterization  used  in  the  Community  Climate 
Model  version  2  (CCM2;  Briegleb  1992).  The  Briegleb  (1992)  solver  in  CCM2 
is  also  available  as  a  shortwave  radiative  transfer  option  in  MM5.  The  CAM2 
model  accounts  for  absorption  by  ozone,  water  vapor,  carbon  dioxide,  and  oxygen, 
scattering  by  molecules,  and  absorption  and  scattering  by  cloud-water  drops. 

Ozone  absorption  was  initially  treated  by  Briegleb  (1992)  as  gray  in  eight 
spectral  intervals  from  0.2-0. 7  pm.  The  absorption  coefficient  values  and  solar 
weights  were  derived  from  the  absorption  cross  section  data  of  King  and  Cast 
(1961),  and  solar  data  used  for  the  weights  were  taken  from  Labs  and  Neckel 
(1968).  Briegleb  (1992)  used  a  trial  and  error  method  to  reduce  the  number 
of  spectral  intervals  from  those  of  King  and  Cast  (1961)  to  the  eight  intervals 
between  0.2-0. 7  pm.  Later,  the  spectral  interval  from  0.35-0.70  pm  was  divided 
at  0.64  pm  to  provide  better  agreement  with  observations  (Collins  1998),  for  a 
total  of  nine  ozone-related  bands. 

Absorption  by  oxygen  accounts  for  the  A-  and  B-oxygen  bands  between  0.35- 
0.70  pm.  Oxygen  absorption  is  treated  as  gray  over  the  interval  and  the  absorption 
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is  scaled  by  the  absorber  amount  uq^  defined  as 


U02 


^/f^oQo2 


.p 

9 


(4.9) 


where  p  is  the  path  pressure,  g  is  gravity,  po  is  the  cosine  of  the  solar  zenith  angle, 
and  qo2  is  the  mass  mixing  ratio  of  oxygen. 

Absorption  by  carbon  dioxide  in  two  near  infrared  bands,  2.63-2.86  pm  and 
4.16-4.55  pm,  is  treated  in  much  the  same  way  as  oxygen  using  data  from  Kiehl 
et  al.  (1985).  Note  that  the  4.16-4.55  pm  band  is  actually  subdivided  into  two 
monochromatic  bands  using  the  weighted-sum-of-gray-gases  method  (Hottel  and 
Sarofim  1967;  Modest  1991). 

As  in  the  GFDL-Eta  model,  the  Briegleb  (1992)  approach  makes  use  of  the 
Lacis  and  Hansen  (1974)  treatment  of  water  vapor  absorption  in  the  near  infrared. 
Briegleb  (1992)  uses  seven  of  the  eight  /c- values  from  Lacis  and  Hansen  (1974),  but 
the  weights  are  empirically  adjusted  to  give  good  agreement  to  LBL  calculations  by 
Ramaswamy  and  Freidenreich  (1991)  of  water  vapor  transmission  (B.  P.  Briegleb 
and  W.  D.  Collins  2004,  personal  communication).  Cloud  optical  properties  are 
taken  from  Slingo  (1989),  where  the  properties  are  parameterized  based  on  liquid 
water  path  and  equivalent  radius  of  the  drop-size  distribution. 

Molecular  scattering  optical  depths  are  taken  from  Frohlich  and  Shaw  (1980), 
and  are  applied  at  the  midpoint  of  each  interval  from  0.20-0.35  pm.  For  the  visible 
band  the  optical  depth  was  evaluated  at  0.495  pm,  the  wavelength  that  gave  the 


best  agreement  with  their  benchmark.  The  molecular  scattering  optical  depth  for 
the  first  quadrature  interval  in  the  near  infrared  was  tuned  to  £t  their  benchmark 
calculations,  while  the  midpoint  of  the  near  infrared  spectral  interval  was  used 
to  calculate  the  molecular  scattering  optical  depths  for  the  remaining  quadrature 
intervals. 

The  CAM  version  3  (CAMS)  radiation  code  (Collins  et  ah  2004)  is  an  update 
of  version  2  with  a  modified  treatment  of  absorption  by  water  vapor  in  the  near 
infrared.  In  CAMS  the  quadrature  weights  are  the  same  as  those  used  by  CAM2. 
The  seven  values  are  determined  using  iterative  exponential  sum  fitting  to  fit  LBL 
calculations.  For  the  hrst  three  /^-values,  Collins  et  ah  (2004)  minimize  differences 
between  flux  divergences  calculated  with  their  parameterized  A:-values  and  weights 
and  those  calculated  by  LBL.  For  the  remaining  four  /c-values,  Collins  et  ah 
(2004)  minimize  differences  between  parameterized  and  LBL-calculated  heating 
rate  profiles. 

The  LBL  calculations  are  based  on  the  HITRAN2000  database  and  its  2001 
update  (Rothman  et  ah  2003)  with  water  vapor  continuum  absorption  modeled 
using  the  CKD-2.4.1  approximation  (Clough  et  al.  1989).  Solar  spectral  informa¬ 
tion  is  taken  from  Kurucz  (1994),  and  Rayleigh  scattering  in  the  near  infrared  is 
calculated  following  Lenoble  (1993).  The  near  infrared  spectral  interval  in  CAM3 
also  includes  the  effects  of  water  vapor  line  absorption  out  to  10  /rm,  although  the 
overall  solar  weight  given  to  the  band  is  unchanged  from  the  CAM2  weight  for 
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the  band  from  0. 7-5.0  /rm.  In  application,  the  only  differences  between  the  CAM2 
and  CAMS  shortwave  radiation  codes  are  the  values  of  the  seven  /c-values  in  the 
near  infrared  and  the  value  of  an  empirical  constant,  6,  used  to  scale  water  vapor 
absorption  to  agree  with  LBL  heating  rates  in  the  stratosphere. 

We  calculate  the  gas  and  cloud  optical  properties  according  to  the  methods 
used  in  the  CAM2  and  CAMS  models  and  pass  them  to  the  Zdunkowski  et  ah 
(1982)  two-stream  solver  to  calculate  heating  rates  for  the  comparisons  in  Section 
4.S. 

4.2.4  Fu  and  Lion  Scheme 

Fu  and  Lion  (1992)  describe  a  correlated  /^-distribution  model  that  divides  the 
solar  spectrum  into  hve  bands  associated  with  water  vapor,  with  a  sixth  band  for 
ozone  absorption  added  later.  The  Fu  and  Lion  (1992)  scheme  uses  a  total  of  44 
(later  54)  quadrature  points.  The  Fu  and  Lion  (1992)  model  bears  a  similarity  to 
the  FSCK  model  for  shortwave  radiation  described  in  Section  2.1.2  because  they 
account  for  the  spectral  dependence  of  the  solar  flux  by  dehning  the  cumulative 
/^-distribution  as 

g{k,p,T)  =  '^Sjgj{k,p,T)Auj/STj,  (4.10) 

j 

where  Sjy  =  Ylj  iSjAuj.  Fu  and  Lion  (1992)  denote  each  of  the  six  spectral  bands 
as  Az/,  which  are  themselves  divided  into  spectral  subintervals  of  width  Auj  taken 
to  be  50  cm“^.  Over  each  of  the  Auj  subintervals  the  solar  flux  is  assumed  to 
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be  constant  with  a  spectrally  integrated  value  of  Sj.  A  cumulative  /^-distribution, 
gj{k,p,T),  is  created  for  each  of  the  subintervals  from  which  the  total  cumulative 
/^-distribution,  g{k,p,T),  for  the  band  is  constructed.  Cumulative  /^-distributions 
are  created  for  11  pressures,  p,  and  3  temperatures,  T. 

The  version  of  the  Fu  and  Liou  (1992)  model  used  for  the  comparison  in  Section 
4.3  is  essentially  that  described  in  the  original  paper,  together  with  the  added  spec¬ 
tral  band  for  ozone  absorption.  We  obtained  the  code  from  the  Spherical  Harmonic 
Discrete  Ordinate  Method  (SHDOM;  Evans  1998)  website.  As  a  consequence,  our 
version  of  the  Fu  and  Liou  (1992)  code  is  not  the  more  recent  version  used  by 
the  Clouds  and  the  Earth’s  Radiant  Energy  System  (CERES)  program  at  NASA 
Langley.  Dehciencies  in  our  version  of  the  Fu  and  Liou  (1992)  parameterization 
include  use  of  a  single  ozone  band  for  the  ultraviolet  and  visible  spectral  regions, 
the  outdated  water  vapor  continuum  absorption  (Roberts  et  ah  1976),  and  missing 
trace  gas  absorbing  species  (Evans  2004,  personal  communication).  We  couple  the 
gas  absorption  k-values  of  Fu  and  Liou  (1992)  with  Rayleigh  scattering  calculated 
according  to  Teillet  (1990)  and  cloud  optical  properties  calculated  according  to 
the  “Mie”  approach  of  Toon  and  Ackerman  (1981)  into  the  two-stream  solver  of 
Zdunkowski  et  ah  (1982)  for  the  comparisons  in  Section  4.3. 
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4.2.5  Kato  Scheme 


The  correlated  /^-distribution  model  developed  by  Kato  et  al.  (1999)  is  a  traditional 
correlated  /^-distribution  method  as  described  in  Section  2.1.1.  The  model  accounts 
for  gas  absorption  by  ozone,  oxygen,  carbon  dioxide,  and  water  vapor,  including 
water  vapor  continuum  absorption,  using  154  quadrature  points  across  32  spectral 
intervals  over  each  of  which  cloud  properties  are  approximately  constant.  Tables  of 
quadrature  values,  similar  to  those  described  in  Chapter  3,  are  generated  for  water 
vapor  and  oxygen,  while  ozone  is  treated  as  a  gray  gas  over  several  of  the  spectral 
intervals  and  carbon  dioxide  and  ozone  are  treated  as  minor  constituents  with 
hxed  concentrations  in  the  /^-distributions.  Using  these  gas  absorption  properties, 
along  with  Rayleigh  scattering  calculated  according  to  Teillet  (1990)  and  cloud 
optical  properties  calculated  according  to  Toon  and  Ackerman  (1981),  we  again 
calculate  heating  rates  with  the  two-stream  solver  of  Zdunkowski  et  al.  (1982)  for 
the  comparisons  in  Section  4.3. 

4.2.6  RRTM  SW  Scheme 

The  shortwave  Rapid  Radiative  Transfer  Model  (RRTM  SW;  Mlawer  and  Clough 
1997;  Mlawer  and  Clough  1998)  is  a  correlated  /;-distribution  model  that  uses  14 
narrow  bands  for  the  shortwave,  with  16  quadrature  points  per  band,  for  a  total  of 
224  radiative  transfer  calculations.  Gas  absorption  values  for  water  vapor,  carbon 
dioxide,  ozone,  methane,  and  oxygen  are  stored  at  59  reference  pressures.  The 
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spectral  bands  “were  chosen  as  a  consequence  of  the  absorption  band  structure  of 
the  absorbing  species  (Mlawer  and  Clough  1997).” 

Applying  an  algorithm  similar  to  the  FSCK  approach,  RRTM  keeps  track  of 
the  “mapping”  of  absorption  coefficients  from  spectral  space  to  g'-space  for  the 
highest  atmospheric  level  in  the  McClatchey  et  al.  (1972)  midlatitude  summer 
atmosphere  as  a  reference  level.  That  mapping  is  then  applied  to  the  Planck, 
or  solar  source,  function  (e.g.,  Kurucz  1994)  and  these  mapped  solar  irradiance 
values  are  summed  and  used  as  weights  for  each  g'-value. 

At  the  reference  level,  this  method  would  be  identical  to  that  described  in 
Section  2.1.2  if  the  mapped  solar  source  function  were  used  instead  of  the  standard 
/^-distribution.  In  order  to  ensure  equivalence  with  the  FSCK  method,  it  would 
be  necessary  to  keep  track  of  the  monochromatic  spectral  intervals  associated  with 
each  mapped  value  of  the  solar  source  function.  This  would  make  it  possible  to 
normalize  the  resulting  energy  distribution  over  the  spectral  interval  to  which  it 
applies.  As  such,  this  method  would  be  exact  for  a  homogeneous  atmosphere 
if  extended  to  the  full  spectrum.  However,  because  the  mapping  of  the  solar 
irradiances  is  only  done  at  the  reference  level,  if  this  technique  were  extended  to 
the  full  spectrum,  the  mapping  would  not  be  correct  at  levels  where  the  sorting 
changes  signihcantly,  i.e.,  where  correlation  breaks  down.  Because  Mlawer  and 
Clough  (1997)  have  implemented  bands  that  take  into  account  correlations  across 
the  shortwave  spectrum,  they  avoid  a  breakdown  of  the  mapping  but  at  a  cost  of 
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many  radiative  transfer  calculations. 

The  RRTM  SW  treatment  of  absorption  by  gases  used  for  this  comparison  is 
the  version  implemented  in  SHDOM  (Evans  1998)  and  differs  from  the  original 
algorithm  of  Mlawer  and  Clough  (1997)  and  Mlawer  and  Clough  (1998).  In  the 
SHDOM  version,  although  the  ^-distribution  mapping  of  the  solar  source  func¬ 
tion  is  implemented,  the  mapped  solar  irradiances  are  summed  for  each  complete 
spectral  band  rather  than  just  for  the  weights  associated  with  each  quadrature 
point.  As  a  result,  some  of  the  correlation  between  solar  flux  and  absorption  is 
lost  (S.  A.  McFarlane  2004,  personal  communication).  Mlawer  and  Clough  (1998) 
show  that  neglecting  the  mapping  results  in  a  substantial  difference  in  energy  ab¬ 
sorbed  by  the  atmosphere.  While  the  mapping  is  not  completely  neglected  for  the 
solar  irradiances  in  this  implementation,  these  differences  in  implementation  and 
their  consequences  must  be  noted  when  interpreting  results  in  the  comparisons  in 
Section  4.3. 

As  for  Fu  and  Liou  (1992)  and  Kato  et  al.  (1999),  we  used  the  gas  absorption 
properties  from  the  RRTM  SW  algorithm,  along  with  cloud  optical  properties 
calculated  according  to  Toon  and  Ackerman  (1981),  to  compute  heating  rates 
with  the  two-stream  solver  of  Zdunkowski  et  al.  (1982).  Furthermore,  nongray 
Rayleigh  scattering  is  normally  combined  with  the  absorption  coefficients  in  the 
RRTM  scheme  (S.  A.  McFarlane  2004,  personal  communication).  The  precise 
methodology  for  combining  Rayleigh  scattering  with  gas  absorption  coefficients 
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and  later  separating  the  two  quantities  for  use  by  the  solver  is  not  clear  from  the 
presentation  of  the  method  (Mlawer  and  Clough  1997;  Mlawer  and  Clough  1998). 
Consequently,  Rayleigh  scattering  for  the  RRTM  scheme  is  gray  in  each  spectral 
band  for  our  comparison,  calculated  according  to  Teillet  (1990).  Mlawer  and 
Clough  (1998)  remark  that  the  use  of  spectrally  averaged  values  of  key  physical 
properties,  such  as  the  Rayleigh  scattering  used  for  this  comparison,  will  result  in 
signihcant  differences  from  the  complete  RRTM  SW  scheme. 

4.3  Comparison  to  State-of-the-Art  Operational 
and  Research  Models 

Figure  4.5  shows  heating  rates  for  the  ICRCCM-III  CLEAR  case  calculated  by 
each  of  the  models  described  in  Section  4.2,  where  height  is  the  vertical  coordinate 
in  the  figure  rather  than  pressure  as  in  the  ICRCCM-III  comparison.  Figure  4.5a 
extends  into  the  mesosphere,  a  part  of  the  atmosphere  not  normally  considered  in 
numerical  weather  prediction.  However,  ozone  must  be  included  as  an  absorbing 
gas  in  the  shortwave,  and  a  model  must  account  for  the  reduction  of  incoming 
solar  radiation  by  ozone  in  the  stratosphere  before  the  radiation  reaches  a  typical 
“model  top”  level  used  for  other  physical  and  dynamical  calculations.  Otherwise, 
too  much  absorption  and  heating  will  occur  in  the  top  layers  of  the  model.  Figure 
4.5a  shows  a  striking  range  of  heating  rate  differences  between  models,  primarily 
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Heating  Rate  (K/day)  %  Error 

Figure  4.5.  Broadband  heating  rates  (K/day)  calculated  by  LBL  (black  solid),  the 
FSCK  approach  (black  dotted),  the  Dudhia  scheme  (light  gray  solid),  the  GFDL-Eta 
scheme  (light  gray  dash-dot),  the  CAMS  scheme  (light  gray  dashed),  the  Fu  and  Lion 
scheme  (dark  gray  solid),  the  Kato  scheme  (dark  gray  dashed),  and  the  RRTM  SW 
scheme  (dark  gray  dash-dot)  for  the  CLEAR  case  of  Barker  et  al.  (2003). 

due  to  differences  in  the  treatment  of  absorption  by  ozone.  Also  note  that  the 
heating  rate  percent  errors  for  all  of  the  models  in  the  comparison  become  large 
and  negative  between  ~12-18  km  as  the  magnitude  of  the  LBL  heating  rate  prohle 
becomes  small,  i.e.,  less  than  0.5  K/day. 
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The  Dudhia  model,  which  does  not  include  absorption  by  ozone,  was  not  de¬ 
signed  to  model  the  atmosphere  above  the  troposphere.  Therefore,  this  model  has 
essentially  100%  heating  rate  errors  relative  to  LBL  calculations  above  the  lower 
stratosphere  (Fig.  4.5a).  Near  the  troposphere  (Fig.  4.5b)  the  Dudhia  model  over¬ 
estimates  heating  by  over  30%  and  underestimates  heating  by  approximately  50% 
near  the  ground.  The  heating  rate  RMS  error  for  the  Dudhia  model  relative  to 
the  LBL  benchmark  over  all  atmospheric  layers  is  6.28  K/day.  Note  that  although 
the  heating  rate  prohles  are  shown  for  the  Dudhia  model  in  the  left  panels  of  Fig. 
4.5a-b  and  Dudhia  model  heating  rate  errors  are  plotted  in  the  right  panel  of  Fig. 
4.5a,  the  Dudhia  model  heating  rate  errors  are  not  plotted  in  the  right  panel  of 
Fig.  4.5b. 

Figure  4.6  shows  the  percent  errors  in  downward  flux  absorbed  at  the  sur¬ 
face  calculated  by  each  of  the  models  relative  to  the  LBL  calculated  value  of 
483.72  Wm-2  for  the  ICRCCM-III  CLEAR  case.  In  numerical  weather  prediction, 
downward  flux  at  the  surface  is  extremely  important  because  it  is  a  major  con¬ 
tribution  to  the  surface  energy  budget.  The  surface  energy  budget  controls,  to  a 
great  extent,  the  operation  of  the  planetary  boundary  layer  model,  which,  in  turn, 
has  a  major  influence  on  the  overall  model  physical  and  dynamical  solution. 

The  Dudhia  model  underestimates  the  downward  flux  at  the  surface  by  nearly 
6%.  While  the  Dudhia  model  errors  in  heating  rate  and  downward  flux  at  the 
surface  are  quite  large  relative  to  the  other  models  in  this  comparison,  the  heating 
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Figure  4.6.  Percent  error  in  model  downward  flux  absorbed  at  the  surface  relative  to 
LBL  calculations  (483.72  Wm“^). 

rates  in  the  troposphere  are  correct  to  within  a  factor  of  2  or  less,  which  may  be 
acceptable  for  some  applications  where  efficiency  has  a  signihcantly  higher  priority 
than  accuracy.  Figure  4.7  shows  that  the  Dudhia  model  requires  ~8  seconds  on  a 
single  1.7  GHz  processor  for  100,000  solutions  of  the  ICRCCM-III  CLEAR  case, 
which  includes  generating  the  optical  properties  and  running  the  radiative  transfer 
solver.  This  is  more  than  four  times  faster  than  the  next  fastest  model. 

The  GFDL-Eta  model  performs  exceptionally  well  in  the  upper  atmosphere 
compared  to  the  other  models,  with  heating  rate  errors  generally  less  than  20% 
above  ~45  km  (Eig.  4.5a).  In  the  lower  atmosphere  the  GFDL-Eta  heating  rates 
are  less  accurate  than  most  of  the  other  models,  underestimating  heating  by  as 
much  as  10-15%  (Fig.  4.5b).  The  heating  rate  RMS  error  for  the  GFDL-Eta 
model  relative  to  the  LBL  benchmark  over  all  atmospheric  layers  is  0.78  K/day. 
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Figure  4.7.  Processing  time  in  seconds  on  a  single  1.7  GHz  processor  for  100,000 
iterative  model  calculations  without  inputs  and  outputs  to  hard  disk.  The  time  spent 
on  creating  the  optical  properties  is  represented  by  the  lightly  shaded  region,  while  the 
time  spent  on  the  solver  is  the  darkly  shaded  region. 


This  deficiency  is  also  evident  in  Fig.  4.6,  which  shows  that  the  GFDL-Eta  model 
overestimated  the  downward  fiux  at  the  surface  by  6.7%  .  Figure  4.7  shows  that 
the  GFDL-Eta  model  was  the  second  fastest  in  this  comparison,  requiring  about 
34  seconds  for  100,000  solutions  of  the  GLEAR  case. 

The  GAMS  model  is  significantly  more  accurate  in  the  troposphere  than  the 
GAM2  model  (not  shown)  for  this  GLEAR  case  due  to  the  use  of  recent  water 
vapor  line  and  continuum  absorption  data.  There  is  no  difference  in  the  treatment 
of  ozone  between  the  GAM2  and  GAMS,  however,  with  upper  atmospheric  heating 
rate  errors  as  large  as  40%  (Fig.  4.5a).  GAMS  heating  rate  errors  in  the  lower 
troposphere  are  small  compared  to  most  other  models  in  the  comparison,  generally 
less  than  5%  below  6  km  (Fig.  4.5b).  The  heating  rate  RMS  error  for  the  GAMS 
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model  relative  to  the  LBL  benchmark  over  all  atmospheric  layers  is  2.00  K/day. 
While  the  improvement  of  CAMS  over  CAM2  in  the  troposphere  is  in  large  part 
due  to  the  use  of  more  recent  information,  the  small  heating  rate  errors  relative 
to  the  LBL  calculations  in  the  comparison  are  an  indication  of  the  care  taken  in 
the  empirical  reevaluation  of  water  vapor  absorption  for  CAMS.  The  CAMS  model 
improvement  over  CAM2  is  evident  in  Fig.  4.6,  which  shows  that  the  CAM2  model 
overestimated  the  CLEAR  case  downward  flux  at  the  surface  by  2.8%  while  the 
CAMS  model  overestimates  it  by  only  O.S%.  Timing  for  the  CAM2  and  CAMS 
models  is  identical,  since  only  eight  parameters  changed  between  the  two  versions. 
The  two  models  each  take  285  seconds  for  100,000  solutions  of  the  CLEAR  case 
(Fig.  4.7),  with  55%  of  that  time  devoted  to  generating  the  optical  properties  and 
the  remaining  45%  of  the  time  required  for  the  radiative  transfer  solver. 

The  version  of  the  Fu  and  Lion  (1992)  model  that  we  have  tested  here  has 
heating  rate  errors  of  nearly  60%  in  the  upper  atmosphere  (Fig.  4.5a)  and  nearly 
1S%  in  the  lower  troposphere  (Fig.  4.5b).  The  heating  rate  RMS  error  for  the 
Fu  and  Lion  model  relative  to  the  LBL  benchmark  over  all  atmospheric  layers  is 
2.8S  K/day.  Because  we  are  not  using  the  most  current  version  of  the  Fu  and 
Lion  (1992)  model,  these  heating  rate  errors  result,  in  part,  from  the  omission  of 
water  vapor  continuum  absorption  and  old  spectral  line  absorption  information 
for  both  water  vapor  and  ozone  absorption.  But,  we  cannot  separate  this  source 
of  error  from  the  choice  of  bands  and  the  number  of  quadrature  points  per  band. 
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The  Fu  and  Lion  (1992)  model  also  overestimates  downward  flux  at  the  surface 
by  3.2%,  or  15.6  Wm“^  (Fig-  4.6)  and  requires  702  seconds  per  100,000  solutions 
of  the  CLEAR  case  (Fig.  4.7).  Of  the  total  processing  time,  44%  was  expended  to 
calculate  the  optical  properties  and  56%  was  used  by  the  solver.  The  generation  of 
optical  properties  requires  5.7  seconds  per  quadrature  point  per  100,000  solutions, 
resulting  mainly  from  the  2-D  interpolation  of  A:-values  stored  on  the  pressure  and 
temperature  grid. 

The  Kato  et  ah  (1999)  model  tested  here  is  based  on  the  same  HITRAN 
database  as  the  LBL  benchmark.  This  fact,  coupled  with  the  large  number  of 
spectral  bands  in  the  Kato  scheme,  leads  to  heating  rate  errors  that  are  less  than 
10%  at  the  top  of  the  atmosphere  (Fig.  4.5a)  and  remain  less  than  30%  through 
the  entire  atmospheric  profile.  Errors  from  the  surface  to  ~6  km  are  less  than  5% 
and  remain  less  than  10%  up  to  10  km  (Fig.  4.5b).  The  heating  rate  RMS  error 
for  the  Kato  model  relative  to  the  LBL  benchmark  over  all  atmospheric  layers 
is  0.85  K/day.  Downward  flux  at  the  surface  is  overestimated  by  3.2  Wm“^  or 
0.7%  (Fig.  4.6).  This  accuracy  results  from  the  use  of  154  quadrature  points 
covering  32  bands,  which  comes  at  the  cost  of  efficiency.  The  Kato  model  required 
3175  seconds  for  100,000  solutions  (Fig.  4.7),  of  which  65%  was  used  to  generate 
the  optical  properties  and  35%  was  needed  to  run  the  solver.  The  generation  of 
optical  properties  requires  13.5  seconds  per  quadrature  point  per  100,000  solutions, 
resulting  mainly  from  the  3-D  interpolation  from  /c-values  stored  on  a  pressure. 
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temperature,  and  water  vapor  mixing  ratio  grid. 

The  RRTM  shortwave  model  heating  rate  errors  are  nearly  50%  in  the  upper 
atmosphere  (Fig.  4.5a),  and  the  heating  rate  prohle  above  ~40  km  bears  a  strong 
similarity  to  that  of  the  Fu  and  Lion  (1992)  model.  In  the  lower  atmosphere, 
however,  heating  rate  errors  are  near  0%  below  4  km,  and  do  not  exceed  10%  in 
the  lowest  12  km  (Fig.  4.5b).  The  heating  rate  RMS  error  for  the  RRTM  model 
relative  to  the  LBL  benchmark  over  all  atmospheric  layers  is  2.28  K/day.  The 
RRTM  shortwave  model  overestimates  downward  flux  at  the  surface  by  0.9%  (Fig. 
4.6)  and  requires  1758  seconds  for  100,000  solutions  (Fig.  4.7).  Of  the  total  time 
used  by  RRTM,  only  9%  was  used  to  generate  the  optical  properties,  while  the 
remaining  81%  was  used  by  the  solver. 

The  time  used  by  the  Zdunkowski  et  ah  (1982)  two-stream  solver  is  directly 
proportional  to  the  number  of  quadrature  points  per  layer.  The  large  proportion 
of  time  used  by  the  solver  with  the  RRTM  shortwave  model-generated  optical 
properties  results  from  the  use  of  224  quadrature  points.  The  amount  of  time 
required  to  generate  the  optical  properties  for  a  correlated  ^-distribution  model 
appears  to  be  related  not  only  to  the  number  of  points,  but  also  the  dimensionality 
of  the  table  in  which  /c-values  are  stored.  In  RRTM  the  A:- values  are  stored  on  a  1-D 
pressure  grid  and  generation  of  the  optical  properties  requires  only  0.7  seconds  per 
quadrature  point  per  100,000  solutions.  To  all  appearances  the  code  to  generate 
optical  properties  for  the  RRTM  shortwave  model  has  been  carefully  optimized. 
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Note  that  we  did  not  use  any  compiler  optimization  options  for  any  of  the  models 
used  in  this  comparison. 

Heating  rates  calculated  by  the  two-band  FSCK  model  remain  within  20%  of 
the  LBL  benchmark  through  the  entire  atmospheric  prohle  (Fig.  4.5a),  with  the 
exception  from  12-18  km  where  LBL  heating  rate  values  are  less  than  0.5  K/day. 
Heating  rate  errors  below  10  km  do  not  exceed  7%  (Fig.  4.5b).  The  heating  rate 
RMS  error  for  the  FSCK  model  relative  to  the  LBL  benchmark  over  all  atmospheric 
layers  is  0.57  K/day.  The  FSCK  model-calculated  downward  surface  fluxes  are 
overestimated  by  5.7  Wm“^  or  1.2%  (Fig.  4.6).  The  FSCK  model  took  391  seconds 
for  100,000  solutions  (Fig.  4.7),  requiring  72%  of  the  time  for  calculating  optical 
properties,  averaging  18.9  seconds  per  quadrature  point  per  100,000  solutions,  and 
18%  for  the  solver. 

The  results  of  the  model  comparisons  for  the  ICRCCM-HI  CLOUD  A  and 
CLOUD  B  cases  are  shown  in  Fig.  4.8.  As  in  Fig.  4.5b,  we  show  heating  rate 
prohles  for  the  Dudhia  model,  but  we  omit  the  Dudhia  model  heating  rate  errors. 
As  mentioned  above,  heating  rates  for  the  GFDL-Eta  model  are  not  shown  for  the 
CLOUD  A  and  CLOUD  B  cases. 

The  Dudhia  model  overestimates  in-cloud  heating  by  nearly  4  K/day  for  the 
CLOUD  A  case  and  by  over  8  K/day  for  the  CLOUD  B  case,  representing  errors 
of  about  66%  and  130%,  respectively.  Again,  these  results  are  due  to  the  use  of 
a  very  simple  empirical  algorithm  that  is  fast  but  not  accurate.  There  is  much 
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Figure  4.8.  Broadband  heating  rates  (K/day)  calculated  by  LBL  (black  solid),  the 
FSCK  approach  (black  dotted) ,  the  Dudhia  scheme  (light  gray  solid) ,  the  CAMS  scheme 
(light  gray  dashed),  the  Fu  and  Lion  scheme  (dark  gray  solid),  the  Kato  scheme  (dark 
gray  dashed),  and  the  RRTM  SW  scheme  (dark  gray  dash-dot)  for  (a)  the  CLOUD  A 
case  and  (b)  the  CLOUD  B  case  of  Barker  et  ah  (2003). 
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less  variability  of  in-cloud  heating  rate  errors  among  the  other  models  in  this 
comparison.  The  CAMS  model,  which  uses  the  Slingo  (1989)  parameterization, 
has  in-cloud  heating  rate  errors  of  about  10%  for  the  CLOUD  A  case  (Fig.  4.8a) 
and  less  than  1%  errors  for  the  CLOUD  B  case  (Fig.  4.8b).  The  CAMS  in-cloud 
heating  rates  are  the  same  as  for  CAM2  (not  shown). 

We  use  500  hPa  transmission-weighted  band  average  “Mie”  theory  cloud  op¬ 
tical  properties  calculated  according  to  Toon  and  Ackerman  (1981)  for  all  three 
correlated  /^-distribution  models,  i.e.,  the  Fu  and  Lion  (1992)  model,  the  Kato 
et  ah  (1999)  model,  and  the  RRTM  shortwave  model  of  Mlawer  and  Clough  (1997) 
and  Mlawer  and  Clough  (1998).  For  the  averaging  of  the  cloud  optical  properties 
we  used  the  spectral  bands  appropriate  for  each  model.  We  also  used  the  Toon 
and  Ackerman  (1981)  method  to  calculate  cloud  optical  properties  for  the  LBL 
benchmark  and  the  nongray  cloud  absorption  for  the  FSCK  model,  although  at  a 
spectral  resolution  of  0.01  /xm.  For  the  LBL  benchmark  and  the  FSCK  model,  we 
then  interpolated  the  0.01  /rm  resolution  cloud  absorption  coefficients  to  the  higher 
resolution  gas  absorption  coefficients  before  using  them.  We  performed  a  similar 
interpolation  of  the  cloud  scattering  coefficients  in  the  LBL  benchmark  calcula¬ 
tions.  As  explained  in  Chapter  2  the  two-band  FSCK  cloud  scattering  properties 
were  tuned  to  minimize  flux  errors. 

As  Fig.  4.8  illustrates,  the  nongray  treatment  of  cloud  absorption  in  the  FSCK 
model  is  an  accurate  alternative  to  treating  cloud  optical  properties  as  gray  in 
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several  narrow  spectral  bands.  The  process  of  averaging  cloud  absorption  over 
typical  narrow  bands,  where  the  absorption  is  rapidly  changing,  clearly  may  intro¬ 
duce  more  error  than  using  directly  nearly  monochromatic  LBL  and  FSCK  cloud 
absorption  coefficients  in  the  calculations. 

Figures  4.5,  4.6,  4.7,  and  4.8  clearly  show  that  the  two-band  FSCK  model  is 
more  accurate  than  most  of  the  empirically  based  models.  Moreover,  the  two-band 
FSCK  approach  is  more  efficient  than  all  of  the  correlated  /^-distribution  models 
used  in  this  comparison  while  rivaling  them  in  accuracy.  The  new  CAMS  model  is 
also  relatively  accurate  and  efficient,  at  least  in  the  lower  troposphere  for  the  test 
cases  in  our  current  set  of  comparisons. 

As  we  showed  in  Chapter  2,  the  FSCK  method  itself  is  not  the  source  of  the 
errors  that  we  see  in  Fig.  4.5.  In  fact,  none  of  the  models  included  in  our  compar¬ 
ison  rival  the  accuracy  of  the  15-quadrature  point  “hardwired”  solution  through 
most  of  the  atmosphere,  as  shown  in  Fig.  4.9.  The  heating  rate  RMS  error  for  the 
15-quadrature  point  “hardwired”  model  relative  to  the  LBL  benchmark  over  all 
atmospheric  layers  is  0.30  K/day.  The  8192-point  “hardwired”  solution,  an  even 
closer  representation  of  the  performance  of  the  FSCK  algorithm  itself,  has  heating 
rate  errors  that  do  not  exceed  5%  below  9  km  or  above  17  km.  The  heating  rate 
RMS  error  for  the  8192-quadrature  point  “hardwired”  model  relative  to  the  LBL 
benchmark  over  all  atmospheric  layers  is  0.16  K/day. 

The  errors  in  our  implementation  of  the  two-band  FSCK  approach  are  intro- 
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Figure  4.9.  Broadband  heating  rates  (K/day)  for  LBL  calculations  (black  solid),  the 
FSCK  8192-point  hardwired  solution  (dark  gray  solid),  the  FSCK  15-point  hardwired 
solution  (light  gray  solid),  the  FSCK  solution  from  the  current  implementation  (black 
dotted),  and  the  CAMS  scheme  solution  (light  gray  dashed)  for  the  CLEAR  case  of 
Barker  et  al.  (2003). 
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duced  in  several  of  the  steps  in  our  model  development  process.  First,  relatively 
small  errors  are  introduced  by  our  choice  of  15  quadrature  points.  These  errors 
are  indicated  by  the  deviation  of  the  15-quadrature  point  “hardwired”  solution 
from  the  8192-quadrature  point  solution  as  shown  in  Fig.  4.9.  Additional  errors 
are  introduced  by  our  choice  of  pressures,  temperatures,  and  mixing  ratios  used  to 
calculate  stored  absorption  coefficients  and  the  method  we  have  chosen  to  interpo¬ 
late  from  the  stored  values  to  the  atmospheric  prohle  of  interest.  These  errors  are 
evident  when  comparing  the  FSCK  tabulated  solution  to  the  15-quadrature  point 
“hardwired”  solution  in  Fig.  4.9. 

However,  even  with  these  errors,  the  two-band  FSCK  approach  is  more  consis¬ 
tently  accurate  through  the  depth  of  the  prohle  than  the  CAMS  model.  As  Fig. 
4.9  illusrates,  the  CAMS  model  overestimates  heating  by  as  much  as  42%  relative 
to  LBL  from  approximately  25-65  km  versus  a  maximum  FSCK  error  of  17%. 
CAMS  calculations  then  underestimate  heating  by  as  much  as  45%  from  about 
10-25  km,  versus  -25%  for  the  FSCK  approach.  Except  for  a  narrow  range  where 
the  CAMS  results  change  from  overheating  to  underheating,  CAMS  again  overes¬ 
timates  heating  relative  to  LBL  calculations  by  more  than  the  FSCK  calculations 
from  about  5-10  km.  From  approximately  S-5  km  the  FSCK  and  CAMS  methods 
produce  similar  heating  rate  errors  relative  to  LBL.  Only  from  the  lowest  model 
layer  to  about  S  km  are  the  CAMS  heating  rates  more  accurate  than  the  FSCK 
calculations  relative  to  LBL  calculations.  The  reduced  accuracy  of  the  FSCK  ap- 
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proach  for  the  lowest  3  km  likely  results  from  poor  resolution  in  pressure  space 
of  the  tabulated  A:-values.  For  the  lowest  atmospheric  layers,  /^-values  are  linearly 
interpolated  between  those  calculated  for  1050  hPa  and  834  hPa. 

In  terms  of  both  accuracy  and  efficiency  the  current  implementation  of  the  two- 
band  FSCK  approach  is  a  viable  alternative  to  any  of  the  state-of-the-art  models 
included  in  our  comparison.  However,  given  the  superior  accuracy  of  the  two-band 
FSCK  approach  demonstrated  by  the  8192-quadrature  point  and  15-quadrature 
point  “hardwired”  solutions,  it  may  be  worthwhile  to  continue  development  of  a 
more  accnrate,  yet  still  efficient,  method  to  select,  store,  and  retrieve  absorption 
coefficients. 
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Chapter  5 


Discussion 

5.1  Summary  and  Conclusions 

A  fundamental  limitation  of  the  correlated  ^-distribution  method  for  radiative 
transfer  is  the  requirement  that  the  Planck  function  be  constant  over  the  spectral 
interval  for  which  the  A:-distributions  are  constructed.  In  Chapter  2  we  presented 
the  FSCK  method  of  Modest  (2003)  which  removes  the  requirement  of  a  constant 
Planck  function  by  creating  Planck  function-weighted  /^-distributions.  We  simpli- 
hed  the  FSCK  equations  for  the  special  boundary-value  problem  of  solar  radiation 
and  implemented  them  using  the  solar  source  function  of  Kurucz  (1994),  rather 
than  the  Planck  function,  to  weight  the  /:-distributions. 

In  Section  2.2,  we  presented  a  series  of  initial  tests  of  the  FSCK  algorithm. 
For  the  hrst  case,  which  included  only  the  effects  of  gaseous  absorption  with  zero 
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surface  albedo,  we  used  the  FSCK  algorithm  to  model  the  full  shortwave  spectrum 
(0.24-4.60  /im)  as  a  single  spectral  interval.  Using  8192  quadrature  points  to  obtain 
a  baseline  estimate  of  the  accuracy  of  the  FSCK  method,  the  heating  rate  RMS 
error  for  the  full  atmospheric  depth  was  approximately  4.2%,  but  was  nearly  20% 
for  the  lowest  14  km.  This  signihcant  tropospheric  heating  rate  error  motivated  us 
to  divide  the  solar  spectrum  into  two  bands  in  order  to  circumvent  a  breakdown 
of  correlation  with  height  in  clear-sky  absorption.  This  application  of  the  FSCK 
method  to  multiple  spectral  intervals  is  similar  to  the  part-spectrum  approach  of 
Modest  and  Riazzi  (2004),  where  absorption  coefficients  are  divided  into  spectral 
regions  based  on  similar  scattering  properties.  After  dividing  the  spectrum  into 
two  bands,  using  8192  quadrature  points  for  each  band,  the  heating  rate  RMS 
error  was  reduced  to  3.0%  for  the  full  atmospheric  depth  and  3.8%  for  the  lowest 
14  km. 

We  then  reduced  the  number  of  FSCK  quadrature  points  to  5  points  in  band 
1  (0.24-0.68  /im)  and  10  points  in  band  2  (0.68-4.60  /im),  and  again  computed 
the  fluxes  and  heating  rates  for  comparison  with  the  8192  quadrature-point  and 
LBL  results.  The  RMS  error  for  clear-sky  broadband  heating  rates  computed 
using  15  quadrature  points  was  7.3%  for  the  full  atmospheric  depth  and  5.1% 
for  the  lowest  14  km.  The  15  quadrature-point  solution  appeared  to  us  as  the 
best  compromise  between  accuracy  and  computational  efficiency  that  we  could 
obtain.  We  computed  fluxes  and  heating  rates  for  cases  with  nonzero  surface 
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albedo,  molecular  (“Rayleigh”)  scattering  and  low-  and  high-altitude  clouds  in 
three  different  atmospheres. 

Numerical  tests  demonstrated  that  both  gray  and  nongray  treatments  of  ab¬ 
sorption  by  cloud  particles  produced  reasonable  in-cloud  heating  rates  relative  to 
LBL  calculations.  When  clouds  were  treated  as  nongray,  in-cloud  heating  rate 
errors  were  approximately  8-12%,  for  both  single  cloud-layer  and  combined  low- 
and  high-cloud  cases.  Gray  cloud  absorption  produced  less  consistent  results  than 
nongray  cloud  absorption,  with  in-cloud  heating  rate  errors  of  less  than  2%  for  both 
single  cloud-layer  cases,  but  6-14%  for  the  combined  low-  and  high-cloud  cases. 
For  both  treatments  of  cloud  absorption,  scattering  by  cloud  particles  was  treated 
as  gray  in  the  two-spectral  bands,  where  we  tuned  the  gray  scattering  parameters. 
That  is,  we  first  estimated  the  scattering  parameters  using  “Mie”  theory,  then  we 
adjusted  them  to  minimize  flux  errors  above  and  below  the  cloud  layer.  Once  com¬ 
pleted  for  a  given  cloud-drop  distribution  effective  radius  and  spectral  band,  the 
gray  scattering  parameters  are  tabulated  and  become  part  of  the  parameterization. 

The  remaining  errors  in  the  two-band  FSCK  8192  and  15  quadrature-point 
results,  when  using  the  nongray  treatment  of  absorption  by  cloud  particles,  may 
be  the  result  of  a  breakdown  of  correlation  with  height,  as  we  observed  in  the 
full-spectrum  clear-sky  case.  Although  the  heating  rate  errors  with  the  nongray 
treatment  are  typical  of,  or  smaller  than,  the  magnitude  of  in-cloud  heating  rate 
errors  of  other  1-D  models,  we  could  subdivide  the  spectrum  into  more  bands  to 
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address  the  possible  breakdown  of  correlation.  Such  a  modihcation  might  improve 
the  accuracy  of  in-cloud  heating  rates,  but  would  require  additional  quadrature 
points  and  increasing  computational  costs. 

In  Chapter  3  we  advanced  the  FSCK  method  from  a  model  that  can  be  applied 
to  an  atmosphere  with  pre-specihed  pressures,  temperatures,  mixing  ratios,  and 
cloud  properties,  as  in  Chapter  2,  to  an  operational  model  that  calculates  fluxes 
and  heating  rates  given  any  input  atmospheric  conditions.  We  developed  a  set 
of  pressures,  temperatures,  and  mixing  ratios  that  we  used  to  calculate  tabulated 
absorption  coefficients  for  both  the  major  and  minor  gas  constituents.  We  tested 
several  different  methods  for  combining  constituent  absorption  coefficients,  either 
during  the  model  development  process  or  on-the-fly  during  model  runtime. 

Of  the  methods  for  combining  gas  absorption  coefficients  that  we  tested,  the 
method  of  precombining  monochromatic  absorption  coefficients  in  spectral  space 
produced  the  smallest  RMS  errors  in  the  heating  rate  prohles.  This  is  not  surpris¬ 
ing  since  the  precombined  gas  method  does  not  assume  that  component  species  are 
statistically  uncorrelated.  Furthermore,  in  order  to  combine  gases  on-the-fly,  the 
absorption  coefficients  for  each  constituent  are  stored  in  tables  of  quadrature  val¬ 
ues.  The  /j-values  are  retrieved  from  the  tables  and  interpolated  to  the  atmospheric 
state  of  interest  before  being  combined,  and  the  process  of  combining  on-  the-fly 
requires  additional  interpolations.  Each  of  these  processes,  i.e.,  quadrature  prior 
to  combining  constituents  and  interpolations  before  and  during  the  combination 
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procedure,  introduces  errors  not  present  in  the  precombined  absorption  coefficient 
method. 

The  serial  method  of  combining  gas  absorption  coefficients  on-the-fly  is  promis¬ 
ing,  especially  regarding  its  flexibility  for  easily  combining  an  arbitrary  number  of 
gas  constituents.  Since  absorption  coefficients  are  tabulated  for  each  individual 
gas,  whenever  the  individual  gas  absorption  coefficients  are  not  affected  by  self¬ 
broadening,  as  in  the  case  of  ozone  and  carbon  dioxide,  they  only  need  to  be 
calculated  as  a  function  of  pressure  and  temperature,  not  mixing  ratio.  The  ab¬ 
sorption  coefficients  of  these  gases  can  then  simply  be  scaled  by  absorber  amount. 
By  reducing  the  number  of  dimensions  in  the  table  of  stored  ^-values,  the  num¬ 
ber  of  interpolations  required  to  retrieve  them  is  reduced.  However,  in  order  for 
the  method  to  be  effective  operationally,  the  number  of  interpolations  must  be 
reduced  further  and  the  accuracy  of  those  interpolations  improved.  One  possible 
way  of  improving  the  quality  of  interpolations  from  /i:mix(fi'mix)  to  k{g)  at  our  chosen 
quadrature  values  would  be  to  use  cubic  spline  interpolations,  although  such  an 
approach  would  reduce  efficiency. 

We  reexamined  the  treatment  of  absorption  by  cloud  particles  as  either  gray 
or  nongray  in  the  context  of  the  two-band  FSCK  approach.  We  were  not  able 
to  hnd  a  single  two-band  set  of  gray  parameters,  i.e.,  extinction  cross  section, 
single-scattering  albedo,  and  asymmetry  parameter,  that  could  accurately  repro¬ 
duce  both  the  CLOUD  A  and  CLOUD  B  LBL  results.  We  were  able  to  tune  the 
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gray  parameters  to  model  accurately  either  of  the  clouds,  but  parameters  tuned 
for  one  cloud  were  inappropriate  for  the  other.  This  result  implies  that  any  oper¬ 
ational  scheme  would  minimally  have  gray  cloud  parameters  that  are  a  function 
of  both  cloud-particle  size  and  their  vertical  location  within  the  atmosphere.  Such 
a  scheme  is  not  practical  to  implement.  As  a  consequence  of  these  results,  we 
did  not  develop  further  the  use  of  gray  cloud  optical  properties  for  the  two-band 
application  of  the  FSCK  method. 

We  determined  that  precombining  nongray  cloud  absorption  coefficients  with 
gas  absorption  coefficients  was  not  operationally  feasible  because  of  the  large  num¬ 
ber  of  ^-values  that  must  be  stored  in  such  an  approach.  However,  combining 
nongray  cloud  absorption  on-the-fly  with  precombined  gas  absorption  proved  to 
be  accurate  and  efficient,  and  this  was  the  approach  that  we  implemented  for  the 
comparisons  performed  in  Chapter  4. 

Comparing  the  current  implementation  of  the  two-band  FSCK  approach  to 
the  15  quadrature-point  “hardwired”  solution,  we  hnd  that  the  accuracy  of  atmo¬ 
spheric  A:-values  interpolated  from  tabulated  A:-values  improves  if  the  tabulated  k- 
values  are  calculated  for  more  pressures  in  the  lower  troposphere.  For  applications 
in  NWP,  one  might  offset  the  additional  stored  A:-values  in  the  lower  troposphere 
by  storing  fewer  values  in  the  upper  atmosphere.  However,  upper  atmospheric 
radiative  transfer  may  not  be  insignificant  in  NWP  applications,  since  incident 
radiation  must  be  attenuated  accurately  in  a  spectrally  dependent  fashion  to  the 
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“model  top”  in  order  to  obtain  accurate  atmospheric  heating  rates  and  surface 
fluxes.  Regardless,  the  accuracy  of  our  FSCK  15  quadrature-point  solution  could 
be  improved  by  incorporating  additional  pressures  into  the  absorption  coefficient 
tables,  with  the  only  penalty  being  a  table  that  is  larger. 

In  Chapter  4  we  evaluated  the  performance  of  the  two-band  FSCK  model  rela¬ 
tive  to  LBL  calculations  and  to  the  performance  of  other  operational  and  research 
radiative  transfer  parameterizations.  The  evaluation  was  carried  out  in  two  stages. 
First,  we  compared  the  FSCK  calculations  to  the  calculations  of  24  one-dimensional 
radiative  transfer  models  that  were  included  in  the  ICRCCM-III  study  of  Barker 
et  ah  (2003).  Overall,  the  two-band  FSCK  approach  had  smaller  errors  relative 
to  the  LBL  calculations  that  we  used  as  a  benchmark  than  the  median  of  the  1-D 
models  included  in  the  ICRCCM-III  study.  In  the  CLEAR  case  the  FSCK  model 
was  closer  than  the  other  1-D  models  to  the  benchmark  for  TOA  albedo,  surface 
absorptance,  atmospheric  absorptance  and  the  majority  of  the  heating  rate  pro¬ 
file  for  both  bands,  i.e.,  the  ultraviolet- visible  band  and  the  near  infrared,  as  well 
as  for  the  broadband  shortwave.  The  median  of  24  1-D  models  had  a  clear-sky 
mean  bias  error  of  -0.27  K/day  relative  to  LBL  benchmark  heating  rates,  while 
the  operational  implementation  of  the  FSCK  approach  had  a  mean  bias  error  of 
-1-0.04  K/day.  Barker  et  ah  (2003)  attributed  the  tendency  of  the  1-D  models  to 
underabsorb  relative  to  the  ICRCCM-III  benchmarks  to  the  fact  that  many  of  the 
1-D  models  were  developed  with  older  water  vapor  line  and  continuum  absorption 
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data.  The  FSCK  model,  however,  is  based  on  the  same  spectroscopic  data  as  the 
LBL  calculations  that  we  used  as  a  benchmark,  thereby  eliminating  a  source  of 
error  present  in  some  of  the  ICRCCM-III  1-D  models.  Another  likely  source  of 
improved  accuracy  in  the  FSCK  approach  is  that  this  approach  explicitly  accounts 
for  spectral  variations  in  solar  irradiance. 

For  the  CLOUD  A  and  CLOUD  B  cases,  FSCK  broadband  in-cloud  heat¬ 
ing  rates  were  closer  to  the  benchmark  than  the  median  of  the  1-D  models.  For 
CLOUD  A  the  FSCK  error  was  6.1%  versus  16.0%  for  the  1-D  model  median,  and 
for  CLOUD  B  the  FSCK  error  was  -3.5%  versus  15.0%  for  the  1-D  model  median. 
However,  although  the  FSCK-calculated  in-cloud  heating  rate  errors  were  small 
compared  to  the  other  1-D  models,  the  FSCK-calculated  atmospheric  absorptance 
was  too  low  compared  to  the  benchmark  and  the  TOA  albedo  was  too  high,  in¬ 
dicating  less  than  optimal  tuning  of  our  gray  cloud  scattering  properties.  Further 
rehnement  of  the  gray  cloud  scattering  parameters  will  require  additional  testing 
on  a  variety  of  cloud  helds  to  assess  the  performance  of  tabulated  values  under 
different  atmospheric  conditions. 

Second,  we  compared  FSCK  calculations  to  those  of  seven  specihc  state-of-the- 
art  operational  and  research  models.  We  described,  in  detail,  the  methods  used  by 
each  of  the  seven  models  to  generate  optical  properties.  In  the  descriptions  of  the 
models  we  illustrated  how  several  operational  scaled  /^-distribution  and  empirical 
models  have  their  basis  in  algorithms  described  by  Lacis  and  Hansen  (1974).  We 
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embedded  the  gas  and  cloud  optical  properties  from  the  FSCK  model,  the  LBL 
model,  and  most  of  the  six  operational  and  research  models  into  the  two-stream 
solver  of  Zdunkowski  et  ah  (1982),  permitting  us  to  make  a  direct  comparison 
of  treatments  of  optical  properties  and  the  heating  rate  prohles  that  result  from 
them. 

The  model  comparisons  make  clear  that  there  is  a  trade-off  between  accu¬ 
racy  and  efficiency.  The  Dudhia  (1989)  scheme  was  the  fastest  of  the  models 
included  in  our  comparison,  requiring  only  7  seconds  for  100,000  solutions  of  the 
CLEAR  case.  However,  the  Dudhia  scheme  produces  a  downward  flux  at  the 
surface  that  underestimates  the  benchmark  by  nearly  6%  and  atmospheric  heat¬ 
ing  rates  that  underestimate  LBL-calculated  values  by  50-100%  through  most  of 
the  atmosphere.  The  clear-sky  RMS  heating  rate  errors  calculated  over  all  at¬ 
mospheric  layers  relative  to  the  LBL  benchmark  for  the  three  empirically-based 
models  were  0.78-6.28  K/day.  The  Kato  et  ah  (1999)  scheme  calculates  relatively 
accurate  heating  rates  through  the  entire  profile,  but  is  computationally  expensive 
because  the  scheme  uses  154  quadrature  points.  RMS  errors  for  the  three  corre¬ 
lated  /^-distribution  models  were  0.85-2.83  K/day.  Although  several  models  out¬ 
performed  the  others  for  particular  parts  of  the  prohle,  such  as  the  CAMS  scheme 
(Collins  et  ah  2004)  and  the  RRTM  SW  scheme  (Mlawer  and  Clough  1997; 
Mlawer  and  Clough  1998)  which  both  had  heating  rate  errors  of  less  than  5%  be¬ 
low  5  km  and  the  GFDL-Eta  scheme  (Lads  and  Hansen  1974)  which  had  heating 
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rate  errors  less  than  20%  above  45  km,  none  performed  as  consistently  accurate 
through  the  depth  of  the  atmosphere  as  the  two-band  FSCK  approach,  which  had 
an  RMS  heating  rate  error  of  0.57  K/day  over  the  full  atmospheric  prohle. 

Explicitly  accounting  for  the  spectral  variability  of  the  solar  source  function 
likely  improves  the  accuracy  of  the  FSCK  approach  over  other  models,  but  the 
contribution  of  this  effect  is  difficult  to  quantify,  especially  since  a  number  of  the 
models  account  for  spectral  variability  to  a  certain  extent.  Heating  rate  errors  for 
the  FSCK  approach  were  less  than  24%  for  the  full  atmospheric  depth,  and  less 
than  7%  in  the  lowest  10  km.  The  superior  accuracy  of  the  FSCK  approach  also 
extended  to  the  CLOUD  A  and  CLOUD  B  cases,  where  the  FSCK  treatment  of 
absorption  by  cloud  particles  as  nongray  outperformed  the  narrow-band  gray  cloud- 
absorption  treatments  of  the  three  correlated  /^-distribution  schemes.  While  FSCK- 
calculated  in-cloud  heating  rate  errors  were  6.1%  and  -3.5%  for  the  CLOUD  A  and 
CLOUD  B  cases,  respectively,  heating  rate  errors  for  the  correlated  /^-distribution 
schemes  ranged  from  23.4-32.2  %  for  the  CLOUD  A  case  and  13.8-17.8%  for  the 
CLOUD  B  case.  The  FSCK  approach  required  391  seconds  for  100,000  solutions 
of  the  CLEAR  case,  while  the  next  fastest  correlated  /^-distribution  scheme  took 
nearly  twice  as  long. 

Although  the  current  implementation  of  the  two-band  FSCK  approach  per¬ 
formed  well  compared  to  the  ICRCCM-III  1-D  models  and  the  six  state-of-the-art 
operational  and  research  models,  there  are  still  errors  that  may  be  eliminated  by 
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further  refining  the  implementation.  None  of  the  models  included  in  the  study  ri¬ 
val  the  accuracy  of  the  FSCK  8192-point  “hardwired”  atmosphere  solution,  which 
had  an  RMS  heating  rate  error  of  only  0.16  K/day.  The  efficiency  of  the  two-band 
FSCK  approach,  with  its  15  quadrature  points,  rivals  the  efficiency  of  the  CAMS 
scheme.  Consequently,  we  will  continue  to  rehne  the  implementation  of  the  FSCK 
scheme  for  shortwave  atmospheric  radiative  transfer  while  the  work  of  coupling  the 
current  FSCK  scheme  with  an  operational  NWP  model  or  GCM  proceeds. 

5.2  Future  Work 

In  Section  4.3  we  saw  that  the  CAMS  model  had  smaller  heating  rate  errors  from 
the  surface  to  about  3  km  relative  to  the  LBL  benchmark  than  the  FSCK  model. 
We  attributed  the  reduced  accuracy  of  the  FSCK  model  over  this  range  of  altitudes 
to  poor  resolution  in  the  tabulated  A:-values  in  pressure  space  for  the  lower  atmo¬ 
sphere.  Because  the  size  of  the  fc-value  table  is  relatively  small  (approximately 
1  megabyte),  we  can  tabulate  fc- values  for  more  pressures  in  the  lower  atmosphere, 
bearing  in  mind  that  for  each  pressure  we  add  to  the  table  we  are  actually  adding 
/j-values  calculated  for  that  pressure  with  every  combination  of  temperature  and 
gas  mixing  ratio  in  the  table. 

An  alternate  method  for  providing  a  table  of  pre-calculated  A:-values  might  be 
the  use  of  a  table  that  could  adapt  dynamically  to  meteorological  conditions  for  a 
particular  geographic  region.  With  such  a  table,  if  A:-values  near  the  atmospheric 
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state  of  interest  exist  in  the  table,  they  are  used  for  the  solution.  If  no  A:-values 
exist  near  the  desired  atmospheric  state,  new  /c-values  are  computed  on-the-fly 
from  a  secondary  database  of  /^-distributions  and  added  to  the  table.  The  table 
size  would  be  managed  by  removing  from  the  table  those  k-vaAues  which  have  not 
been  used  over  a  specihed  number  of  model  iterations. 

Another  aspect  of  the  FSCK  operational  model  that  may  beneht  from  addi¬ 
tional  research  is  the  treatment  of  cloud  optical  properties.  In  Section  4.1  we  noted 
that  a  likely  cause  of  underestimated  TOA  albedo  by  the  FSCK  model  was  less 
than  optimal  tuning  of  the  gray  cloud  scattering  parameters.  Since  multiple  com¬ 
binations  of  cloud  scattering  cross  sections  and  asymmetry  parameters  can  result 
in  the  same  tuned  solution,  we  must  empirically  determine  the  single  combination 
of  these  parameters  that  is  most  robust  for  a  variety  of  real  cloud  prohles.  Neces¬ 
sary  tests  of  optimally  tuned  parameters  would  include  variations  of  water  content 
and  cloud  physical  thickness  for  both  liquid  and  ice  clouds. 

Inherent  to  all  radiative  transfer  parameterizations  in  NWP  models  is  the  ap¬ 
proximate  representation  of  clouds,  as  well  as  clear  sky,  in  the  models  by  discrete 
vertical  layers.  Radiative  heating  rates  are  strongly  dependent  upon  the  vertical 
thickness  of  model  layers;  that  is,  thin  layers  will  have  signihcantly  higher  heating 
rates  than  thick  layers  with  the  same  optical  thickness.  Future  tests  of  FSCK 
model  treatment  of  clouds  should  examine  how  the  FSCK  model-calculated  heat¬ 
ing  rates  compare  to  those  of  other  models  as  the  vertical  resolution  of  cloud  layers 


131 


changes. 

We  also  saw  in  Section  4.1  that  the  LBL  in-cloud  heating  rates  increased  by 
1.5-2. 5%  when  the  near-infrared  spectral  interval  was  extended  from  4.6  pm  to  5.0 
pm.  This  increase  in  the  heating  rate  is  a  result  of  strong  cloud  absorption  in  the 
spectral  range  4. 6-5.0  pm,  even  though  this  spectral  interval  contains  only  0.1%  of 
the  total  solar  incident  irradiance.  Consequently,  future  development  of  the  FSCK 
method  should  extend  the  solar  spectral  interval  farther  into  the  near-infrared. 

In  addition  to  the  shortwave  version  of  the  FSCK  approach  developed  in  this 
work,  the  full  formulation  of  the  FSCK  method  of  Modest  (2003)  presented  in 
Chapter  2.1.2  describes  a  model  that  can  be  used  for  longwave  atmospheric  ra¬ 
diative  transfer.  Such  a  longwave  model  is  in  the  planning  stage.  Most  likely,  a 
longwave  implementation  of  the  FSCK  method  will  compare  as  favorably  to  op¬ 
erational  and  research  longwave  models  as  the  shortwave  version  has  compared  to 
state-of-the-art  shortwave  models  in  terms  of  accuracy  and  efficiency. 

As  with  the  state-of-the-art  shortwave  models,  there  are  a  variety  of  methods 
used  in  longwave  models.  Fast  techniques  that  rely  heavily  on  empirical  meth¬ 
ods  include,  for  example,  the  Dudhia  longwave  model  (Dudhia  1989)  and  the 
GFDL-Eta  longwave  model  (Fels  and  Schwarzkopf  1975;  Schwarzkopf  and  Fels 
1991).  Implementations  of  the  correlated  /^-distribution  approach  include  that  of 
Fu  and  Lion  (1992),  which  uses  67  quadrature  points  over  12  spectral  intervals,  and 
the  longwave  version  of  RRTM  (Mlawer  et  ah  1997;  Mlawer  and  Clough  1997; 
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Mlawer  and  Clough  1998),  which  uses  240  quadrature  points  over  15  spectral  inter¬ 
vals.  We  anticipate  that  the  FSCK  approach  can  achieve  similar  or  more  accurate 
heating  rates  across  the  longwave  spectrum  with  signihcantly  fewer  calculations 
than  existing  correlated  /^-distribution  methods,  but  this  remains  to  be  demon¬ 
strated. 
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